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ABSTRACT 


A  computer  program  previously  developed  to  estimate  the  motions  and 
accelerations  of  a  planing  craft  in  regular  waves  was  modified  and  extended 
to  compute  the  corresponding  motions  in  random  or  irregular  waves.  Ten 
regular  waves  with  random  phase  were  combined  to  represent  the  random 
seaway.  The  amplitudes  and  frequencies  that  were  selected  represent  the 
energy  distribution  of  a  Pierson-Moskowitz  spectrum  for  a  fully  developed 
sea.  A  comparison  of  computed  results  with  experiments  indicate  that 
the  computer  program  can  predict  craft  behavior  with  reasonable 
quantitative  accuracy  in  moderate  operating  conditions.  In  severe 
operating  conditions,  however,  the  amplitudes  of  the  computed  vertical 
accelerations,  which  include  impacts,  are  one-half  of  the  experimental 
value. 

ADMINISTRATIVE  INFORMATION 

This  work  has  been  authorized  by  the  Naval  Material  Command  (08T2); 
under  Program  Element  625  43  N,  Task  Area  ZF43-421-001,  administered  by 
the  Ship  Performance  Department,  High  Performance  Vehicle  Program  Office, 
Code  1512. 
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INTRODUCTION 


In  a  previous  study*  a  computer  program  was  developed  to  estimate 
the  motions  and  accelerations  of  a  planing  craft  in  regular  waves.  As 
a  logical  extension  of  this  work,  the  program  was  modified  to  compute 
the  motion  of  the  craft  in  random  or  irregular  waves. 

Since  the  mathematical  model  is  nonlinear,  the  computations  are 
made  in  the  time  domain.  Ten  regular  waves  are  combined  with  random  phase 
to  represent  the  random  seaway.  The  amplitudes  and  frequencies  are  ad¬ 
justed  to  conform  to  the  energy  distribution  in  a  Pierson-Moskowitz  fully 
developed  sea. 

The  mathematical  model  was  developed  for  a  V-shaped  prismatic-body 
with  hard  chines  and  constant  deadrise  planing  at  constant  speed.  The 
thrust  and  the  friction  drag  forces  are  assumed  to  act  through  the  center 
of  gravity.  The  vertical  components  of  the  thrust  and  fiction  drag  are 
also  assumed  to  be  negligible  in  comparison  to  the  hydrodynamic  forces. 

The  mathematical  formulation  is  analogous  to  low-aspect- ratio  wing 
theory  with  provisions  for  including  hydrodynamic  impact  loads,  essentially 
a  strip  theory.  Surface  wave  generation  and  forces  associated  with 
unsteady  circulatory  flow  are  neglected,  and  the  flow  is  treated  as 
quasi-steady.  The  mathematical  formulation  is  an  empirical  synthesis  of 
several  theoretically  derived  flows  describing  the  overall  craft  hydro- 
dynami cs . 

MATHEMATICAL  FORMULATION 

EQUATIONS  OF  MOTIONS 


The  aquations  of  motion  for  a  planing  craft  restricted  to  pitch  e, 


heave  z^,  and  surge 


can  be  written  as 
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MXcq  =  Tx  -  N  sin  e  -  0  cos  e 
MzCG  "  t2  "  N  cos  e  +  D  sin  9  +  W 

Ie  =  Nx  -  Dx  .  +  Tx 
cap 


where  M  is  mass  of  craft 

I  is  pitch  moment  of  inertia  of  craft 
N  is  hydrodynamic  normal  force 

0  is  friction  drag 

W  is  weight  of  craft 

Tx  is  thrust  component  in  x  direction 

x  is  distance  from  center  of  gravity  (CG)  to  center  of  pressure 

v 

for  normal  force,  positive  forward 

xd  is  distance  from  CG  to  line  of  action  for  friction  drag  force 

x„  is  moment  arm  of  thrust  about  CG. 

P 

Motions  are  measured  relative  to  a  fixed  coordinate  system  with  the  x  axis 
located  in  the  undisturbed  free  surface  pointing  in  the  direction  of  travel 
and  the  z  axis  pointing  downward. 

Since  the  perturbation  velocities  in  the  forward  direction  are  small 
in  comparison  to  the  speed  of  the  craft,  the  equations  of  motion  can  be 
simplified  by  neglecting  them  and  by  setting  the  forward  velocity  equal 
to  a  constant,  i .e. , 

*CG  =  CONSTANT 

Furthermore,  if  it  is  assumed  that  the  vertical  components  of  the  thrust 
and  friction  drag  forces  are  small  in  comparison  to  the  hydrodynamic  forces 


i 
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and  that  the  total  thrust  and  friction  drag  forces  are  acting  through  the 
center  of  gravity  (so  as  to  produce  no  moments)  the  equations  of  motion 
can  be  written  as 

XCG  =  0 

Mzcg  =  -Ncose  +  W 
Ie  =  Nx 

c 

A  so  called  "strip  theory"  is  used  to  obtain  the  hydrodynamic  force 
acting  on  the  body  by  integrating  the  2-D  hydrodynamic  forces  normal  to  the 
baseline  over  the  wetted  length  of  the  body.  A  body  coordinate  system 

(S,e)  with  its  origin  at  the  CG  and  the  ?  axis  pointing  forward  parallel 
to  the  baseline  of  the  body  as  shown  in  Figure  1  is  used  to  facilitate 
this  integration. 

The  normal  hydrodynamic  force  per  unit  length  f,  acting  at  a  section, 
is  assumed  to  be  proportional  to  the  rate  of  change  of  momentum  associated 
with  an  added  mass  term  and  the  cross  flow  drag,  i.e. 

f  ■  of  K  v>  +  C0,c‘>bv2 

where  V  is  the  velocity  in  plane  of  the  cross  section  normal  to  the 
baseline 

ma  is  the  added  mass  associate  with  the  section  form 

CD,c  ls  the  crossflow  drag  coefficient 

p  is  the  density  of  the  fluid 

b  is  the  half  beam 

Expanding  the  momentum  term  results  is 

^  («,*)  -  ma}  ♦  »Aa  .  (m>v)d| 

where  5  is  the  body  coordinate  parallel  to  the  baseline-  see  Figure  1. 


The  last  term  on  the  right-hand  side  of  the  above  equation  takes  into 
account  the  variation  of  the  section  added  mass  along  the  hull.  This 
contribution  can  be  visualized  by  considering  the  2-D  flow  plane  as  a  sub¬ 
stantive  surface  moving  past  the  body  with  velocity  U  =  -dc/dt  tangent 
to  the  baseline.  As  the  surface  moves  past  the  body,  the  section  geometry 
in  the  moving  surface  may  change  with  a  resultant  change  in  added  mass. 

This  term  exists  evan  in  steady-state  conditions  and  is  the  lift-producing 
factor  in  low-aspect-ration  theory. 

For  a  V-shaped  wedge  the  2-D  added  mass  is  defined  as 
ma  =  **a  f*5  U»t.t) 

where  k  is  the  added  mass  coefficient  (assumed  to  be  1  in  this  study)  and 

a 

b  is  the  wetted  half  beam.  Once  the  chine  becomes  wetted  the  beam  is 

assumed  constant  regardless  of  depth  of  penetration. 

Cross-flow  drag  for  a  V  section  with  separation  at  the  chine  is  assumed 

equal  to  the  drag  of  a  flat  plate  (Cn  _  *  1.0)  corrected  by  the  Bobyleff 

flow  coefficient  approximated  by  cos  e,  i.e., 

Cn  =  1.0  cos  6 
D,C 

The  Bobyleff  flow  coefficient  is  the  theoretical  ratio  of  the  pressure 
on  a  V-section  to  that  experienced  by  a  flat  plate  for  Helmholtz- type 
flow. 

The  same  approximation  is  used  for  estimating  the  drag  coefficient 
for  nonwetted  chine  sections,  using  the  instantaneous  value  of  the  half¬ 
beam  at  the  free  surface. 

An  additional  force  acting  on  the  body  Is  the  buoyancy  force  fg. 

This  force  is  assumed  here  In  to  act  in  the  vertical  direction  and  to 
be  equal  to  the  static  buoyancy  force  multiplied  by  a  correction  factor, 
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i  .e. , 

fB  *  ap9AR 

where  AR  is  the  cross-sectional  area  of  the  section,  and  a  is  a  correction 

factor.  The  full  amount  of  the  static  buoyancy  is  not  realized  because 

at  planing  speeds  the  water  separates  from  the  transom  and  chines,  reducing 

the  pressure  at  these  locations  to  atmospheric  or  less  than  the  equivalent 

hydrostatic  pressure.  A  greater  reduction  is  realized  in  the  buoyancy 

moment  because  of  the  corresponding  shift  in  the  center  of  pressure. 

2 

Shuford  in  his  work  on  steady-state  planing  recommended  a  factor  of 
one- half  to  obtain  the  correct  buoyancy  force.  In  the  following  compu¬ 
tations,  the  buoyancy  force  was  corrected  by  a  factor  of  one-half,  i.e., 
a  =  1/2.  The  buoyancy  moments,  computed  as  the  static  buoyancy  force 
multiplied  by  its  corresponding  moment  arm,  was  corrected  by  an  additional 
factor  of  one-half  to  obtain  the  proper  mean-trim  angles. 

Integrating  the  2-D  hydrodynamic  force  over  the  wetted  length  of  the 


craft  U)  and  taking  the  component  in  the  z  direction  results  in 
-N  cos  0  =  Fz(t)  =  f f  cos  8  d£  +  f fm  it 

*  {/K(WO  +  m,(M)V(*,t) 

“U(U)a*  K<t.t)v(*,t)i 

+  Co>c(*,t)pb({,t)V2(*,t))cos0d* 


Similarly,  the  hydrodynamic  moment  about  the  CG  is  obtained  by  integrating 


the  product  of  the  normal  force  per  unit  length  by  the  corresponding 


moment  arm. 
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F*"i  -.A  cos  mi 

mt«,t)V(*,t)  +  iha(*,t)V({,t) 


+  apA^osfl)  $d| 


Wave  excitation  enters  into  the  above  equations  through  the  geometrical 
properties  of  the  wave,  altering  the  wetted  length  and  draft  of  the 
craft,  and  by  the  vertical  component  of  the  wave  orbital  velocity  at  the 
surface  wz,  altering  the  normal  velocity  V.  Diffraction  has  been  neglected. 

The  horizontal  component  of  orbital  veolcity  is  neglected,  since 
it  is  assumed  small  in  comparison  with  the  forward  speed  The 

velocities  U  and  V  may  then  be  written  as 

U  =  XCG  c°s  8  -  (zcc  -wz)  sin  8 
v  =  xCG  sin  9  -  de  +  (zcc  -wz)  cos  8 

The  depth  of  submergence  h  of  the  body  at  any  point  P(c,c)  may  be  deter¬ 
mined  by 

h  *  zcc  -  $  sin  8  +  f  cos  8  -  r 

where  r  is  the  instataneous  value  of  the  wave  elevation  directly  above 
the  point  measured  positive  downward. 

A  more  detailed  derivation  of  the  above  integrals  for  the  hydrodynamic 
force  and  moment  is  presented  in  reference  (1).  Although  the  hydrodynamic 
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forces  and  moment  require  Integration  over  the  wetted  length,  which 
may  vary  with  time,  the  resulting  equations  of  motion  can  be  integrated 
in  the  time  domain  using  numerical  method  such  as  the  Runge-Kutta  Merson 
integration  routine  used  in  these  studies. 

REPRESENTATION  OF  THE  SEAWAY 

The  seaway  in  general  can  be  represented  by  an  infinite  sum  of  sine 
waves  with  random  phase.  In  these  studies,  for  the  sake  of  computational 
economy,  the  seaway  is  represented  by  the  discrete  sum  of  ten  harmonic 


waves  with  random  phase 


=  £l  ri  '“[kf  (**c, t)*^] 


where  r^  is  the  wave  amplitude 


is  the  wave  number 


c.  is  the  wave  celerity 

and  cr..  is  the  random  phase  angle  of  the  ith  wave  component. 

Note  that  at  point  P(e,c)  on  the  craft 
x  =  xCG  +  ccos  e  +  c  sin  e 

and  XCG  =  ^XCGdt 

Each  frequency  and  wave  amplitude  is  weighted  in  accordance  with  the 
energy  distribution  in  a  Pierson-Moskowitz  spectrum  for  a  fully  developed 
sea.  The  Pierson-Moskowitz  formulation  for  a  continuous  spectrum  can 


be  written  as 


S(w)  *  e 


4flf  0-B/u,4 


where  A  *  8.1  X  10 

g  *  gravitational  constant 

■  ■  ftU* 

The  constant  B  is  also  related  to  the  peak  frequency  Up)  of  the  spectrum 


i 


by 
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which  can  be  confirmed  by  differentiating  the  spectrum  formulation  and 
setting  it  equal  to  zero.  Normalizing  the  frequencies  by  the  peak 
frequency  leads  to  a  nondi mens Iona!  spectrum  S  which  is  related  to  the 
dimensional  spectrum  by  A 

5(n)  « 

where  a  =  u/up 
and  where 

The  discrete  frequencies  representing  the  spectrum  varied  from 

n  =  .80  to  2.6  in  nearly  equal  increments  An=0.2.  A  slight  random 

perturbation  is  given  to  each  frequency  to  avoid  precise  integer  multiple 

frequencies, thereby  increasing  the  fundamental  repetition  period  of  the 

computed  time  history.  Each  discrete  amplitude  is  adjusted  so  that  its 

energy  corresponds  to  that  contained  in  a  band  width  (An)  centered  about 

its  frequency  in  the  continuous  spectrum  i.e., 

+  An/2 

r?  •  ^3  •  /  *h>dB 

n^An/2 

The  band  widths  are  equally  spaced  between  frequencies  except  for  the 
first  and  last  frequencies  which  lump  all  of  the  remaining  energy  at  the 
beginning  and  end  of  the  spectrum.  Table  1  presents  a  list  of  the 
amplitude  for  each  nondi mens ional  frequency  in  terms  of  the  significant 
wave  height. 

COMPARISON  OF  COMPUTED  RESULTS  WITH  EXPERIMENTS 
Computations  of  pitch  and  heave  motions  and  bow  and  CG  vertical 
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^  S(«J*fl«»p)  =  i  e"5/4n 
Hl/3  n 


r * ■  ■'» x ».  t\  - - 1 itnnfltiTW 


accelerations  were  made  using  the  computer  program  (see  Appendix)  for 

3 

comparison  with  the  model  experiments  of  Fridsma  .  Fridsma  tested  a  series 
of  cons tan t-deadrise  models  of  various  lengths  in  irregular  waves  to 
determine  the  effects  of  deadrise,  trim,  load,  speed  and  sea  state  on  the 
added  resistance,  heave  and  pitch  motions  and  vertical  accelerations  at  the 
bow  and  CG.  Figure  2  shows  the  lines  of  the  prismatic  models.  The 
computations  were  made  with  the  Centers'  Control  Data  Corporation  6700 
computer  system.  A  listing  of  the  computer  program  is  presented  in  the 
Appendix. 

Table  2  presents  the  characteristics  of  the  model  craft  for  those 
conditions  selected  for  the  comparison.  The  number  of  computer  runs 
was  kept  to  a  minimum  for  economic  reasons.  Approximately  one  minute 
of  central  processor  time  was  required  for  every  second  of  data  using  model 
scale  dimensions,  and  approximately  40  seconds  of  model  scale  data  was 
required  to  obtain  100  cycles  of  amplitude  data. 

The  output  of  the  program  is  the  time  histories  of  the  pitch  and 
heave  motions  and  the  bow  and  CG  accelerations.  Sample  plots  of  the  out¬ 
puts  are  shown  in  Figures  3  and  4.  Procedures  required  for  processing 
and  analysing  the  data  are  not  a  part  of  this  study.  In  order  to  facili¬ 
tate  the  comparison,  the  analysis  procedures  followed  in  this  report  are 
those  used  by  Fridsma  for  his  experimental  data.  The  amplitudes  (maxima 
or  minima)  of  the  pitch  and  heave  motions  about  the  mean  are  assumed  to 
be  described  by  the  so  called  "Generalized  Rayleiqh  Distribution."  i.e., 

P(y>  ■  7ET  «1/2yV  *  (1  -  e2)1/zy  e-l/2x2dx 

where  y  =  maximum  or  minimum  (absolute  values)  normalized  by  the  standard 
deviation. 
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and  r'  =  ratio  of  negative  maxima  to  total  maxima  or  positive 
minima  to  total  minima. 

To  fit  the  data  to  this  distribution  the  'crests  (maxima)  and  troughs 
(minima)  relative  to  the  mean  value  are  first  determined  from  time  histories 
of  the  motions.  The  mean  value  is  defined  as  halfway  between  the  average 
crest  and  average  trough  value.  The  crests  or  trough  data  (X^ )  are  sorted 
in  ascending  order  and  grouped  into  fifteen  intervals.  At  the  same  time, 
the  proportion  r'  of  negative  maxima  to  total  maxima  or  positive  minima  to 
total  minima  is  determined.  The  cumulative  frequency  and  corresponding 
probability  that  a  crest  or  trough  is  less  than  or  equal  to  the  interval 
value  (X.)  is  then  computed.  From  the  probability  and  r'  values,  the 
theoretical  value  of  the  normalized  amplitude  (y)  is  calculated.  A  plot 
of  versus  corresponding  y  values  is  compared  with  a  line  drawn  through 
x  =  y  =  0,  and  the  point,  x  -  x,  y  -  y  -  /n72(l-2r')  which  is  indicative 
of  the  fit  of  the  theoretical  distribution  function  to  the  data.  The 
values  x  and  y  are  the  observed  average  value  of  the  first  moment  and  the 
theoretical  average  value  (normalized)  respectively. 

Figures  5  and  6  show  typical  examples  of  such  plots  for  the  pitch 
and  heave  crests  and  troughs.  As  can  be  seen  in  these  figures,  the  data 
fit  the  assumed  probability  function  reasonably  well,  but  it  is  also  quite 
possible  that  some  other  distribution  might  fit  the  data  better. 

The  acceleration  data  were  assumed  to  follow  a  simple  exponential 
distribution.  For  this  distribution,  the  probability,  P  of  the  accelera¬ 
tion  peak  n  being  less  than  a  given  value  is 

P(n)  =  1  -  e"n/n 
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where  77  =  average  peak  acceleration. 

Only  the  negative  peak  accelerations  (impact  spikes  as  well  as  wave 
induced)  were  analyzed.  The  data  were  sorted  and  grouped  into  fifteen 
intervals  similar  to  the  motion  data  analysis  and  the  probability  was 
plotted  with  respect  to  n  on  inverted  semilog  paper.  For  a  good  fit,  the 
data  should  follow  a  straight  line  through  the  point  (P  =  0.368,  n  =  n) 
and  the  origin.  Figure  7  shows  a  sample  of  the  acceleration  data  plotted 
in  the  above  manner.  An  exponential  probability  function  appears  to  be  a 
good  fit  to  the  data. 

Table  3  presents  a  comparison  of  the  computed  motions  and  accelerations 
with  the  corresponding  experimental  results.  The  computations  were  made 
for  a  craft  with  a  length  to  beam  ration  (L/b)  of  5,  a  deadrise  angle  (e) 
of  20  degrees  and  a  speed  length  ratio  (V//T)  of  6  for  several  sea  states 
with  significant  wave  height  to  beam  ratios  (H^^b)  of  .222,  .444  and 
.666  which  would  correspond  to  Sea  states  2,  3,  and  4  respectively  for 
a  50  foot  craft.  Computations  were  also  made  in  sea  condition  with 
Hyg/2b  of  .444  for  a  speed  length  ratio  of  4  at  10,  20  and  30  degrees 
deadrise  angles. 

The  tabulated  values  for  the  heave  and  pitch  are  those  with  a  50 
percent  and  those  with  a  90  percent  probability  of  not  being  exceeded. 

Heave  is  nondimensionalized  by  the  beam.  The  values  for  the  bow  and 
CG  accelerations  are  the  average  values  of  the  negative  peaks. 

Other  statistical  variables  such  as  the  1/3  or  1/10  highest  values 
can  be  computed  from  the  specific  probability  distribution.  For  the 
assumed  distribution  of  the  motion  amplitudes,  the  1/10  highest  value  is 
related  to  the  90%  probability  value  by  the  ratio  of  yj/jg/ygo*  wlv*ch  1S 


approximately  1.33  over  the  r'  value  range  measured.  The  1/10  highest 
accelerations  from  the  exponential  distribution  is  3.30  n. 

Plots  of  the  data  in  Table  3  are  also  presented  in  Figure  8  through 
13.  Figures  8  through  10  show  respectively  a  comparison  between  the 
computed  and  experimental  results  of  the  variation  in  the  90  percent 
probability  values  for  the  pitch  and  heave  motions,  and  the  average  values 
for  the  bow  and  CG  accelerations,  with  significant  wave  height  to  beam 
ratio.  Figures  11  through  13  show  similar  plots  for  the  variation  with 
deadrise  angle. 

The  pitch  data  in  the  figures  show  that  while  the  computed  troughs 
(bow  down)  are  in  reasonable  agreement  with  the  experiment,  the  computer 
crests  are  lower  than  the  experiment.  The  heave  exhibits  the  same  trend 
with  the  computed  crests  (CG  down)  being  in  reasonable  agreement  with  the 
expeirment  and  troughs  being  lower  than  the  experiment.  Furthermore,  the 
90  percent  probability  values  for  the  pitch  and  heave  crests  and  troughs 
for  the  computer  model  are  about  equal  in  magnitude;  whereas,  the  experimental 
model  values  in  the  pitch  bow  up  and  heave  CG  up  direction  are  greater. 

It  appears  that  the  experiment  model  exhibits  more  nonlinearity  than  the 
computer  model. 

The  computed  acceleration  data  for  the  bow  and  CG  are  generally  lower 
than  the  corersponding  experiment  data.  Figures  10  and  13  show  that  the 
computed  accelerations  differ  by  15  percent  to  50  percent  of  the  corre¬ 
sponding  experiment  .values  with  the  largest  differences  occuring  in  the 
more  severe  conditions  where  the  accelerations  are  extremely  high. 

For  example,  at  the  ten  degree  deadrise  angle  condition  presented 
in  Table  2,  the  average  1/10  highest  bow  acceleration  corresponds  to  about 


24  g's  for  the  experiment  and  11  g's  for  the  computed  value.  It  is 
doubtful  that  any  operational  boat  would  be  driven  to  such  extreme  conditions. 

The  characteristics  of  the  experiment  model  motions  appear  to  be 
slightly  different  from  those  of  the  computer  model  in  very  extreme 
conditions.  The  experiment  model  experienced  larger  pitch  bow  up  and  heave 
motion  than  the  computer  model,  which  probably  resulted  in  larger  impact 
accelerations.  This  does  not  completely  explain  the  differences  between 
the  experiment  and  computed  vertical  accelerations  which  may  also 
reflect  a  deficiency  in  the  mathematical  representation  of  the  impact 
phenomenon,  and  perhaps  the  seaway  as  well. 

In  less  severe  conditions  the  agreement  between  computed  and  experiment 
results  is  better.  Good  results  can  be  expected  for  speed  length  ratios  up 
to  approximately  6  in  a  seaway  with  significant  wave  height  to  beam  ratio 
of  0.222  (State  2  sea  for  a  50  foot  craft).  For  a  significant  wave  height 
to  beam  ratio  of  0.444  (State  3  for  a  50  foot  craft),  the  calculations 
could  probably  be  used  to  predict  reasonable  quantitative  results  up  to  a 
speed  length  ratio  of  4.  In  more  extreme  conditions  the  comouted  results  are 
less  accurate  quantitatively,  but  are  still  indicative  of  gross  trends. 

SUMMARY  AND  CONCLUSIONS 

A  computer  program  was  developed  to  compute  the  motions  and  accelerations 
of  simple  prismatic  planing  craft  in  head  irregular  waves.  This  was 
achieved  by  incorporating  irregular  waves  into  an  existing  program  for 
computing  the  motions  in  regular  waves.  The  irregular  waves  were  synthe¬ 
sized  by  combining  ten  regular  waves  with  random  phase  and  with  frequencies 
and  amplitudes  weighted  to  represent  a  Pierson-Moskowitz  spectrum  for 
fully  developed  seas. 


Computations  were  made  for  a  craft  with  a  length  to  beam  ratio  of  5, 
a  deadrise  angle  of  20  degrees  and  a  speed  length  ratio  of  6  for  several 
sea  states,  and  in  a  single  sea  state  for  a  speed  length  ratio  of  4  with 
10,  20  and  30  degrees  deadrise  angles.  The  results  were  compared  with  the 
experiment  results  of  Fridsma.  First  the  probability  distributions  of  the 
crests  and  troughs  of  the  motions  and  acceleraions  were  examined  to  deter¬ 
mine  whether  or  not  they  were  the  same  as  those  obtained  by  the  experiments. 
It  was  found  that  the  fit  of  the  pitch  and  heave  crests  and  troughs  to  the 
"Generalized  Rayleigh  Distribution"  which  was  used  in  the  experiment  data 
analysis,  was  acceptable  for  the  computed  data,  but  no  attempt  was  made 
to  fit  the  data  to  other  types  of  distribution  which  might  have  fitted 
better.  The  computed  acceleration  data,  fitted  an  exponential  type  of 
distribution. 

A  comparison  of  the  motions  showed  that  the  computed  pitch  troughs 
were  in  good  agreement  with  the  experiment,  but  the  crests  (bow  up)  were 
lower  than  the  experiment  values.  The  heave  exhibited  the  same  trend,  with 
the  troughs  (CG  heave  up)  for  the  experiment  being  higher  than  the  computed 
values  while  the  crests  were  in  good  agreement. 

Computed  vertical  accelerations  at  the  bow  and  CG  were, for  certain 
conditions,  much!)ower  than  the  experiment  values.  At  some  conditions 
examined,  the  computed  accelerations  were  about  half  of  the  comparable 
experiment  values.  This  occurred  at  the  more  extreme  operating  condition 
where  very  large  acceleations  were  experienced.  For  example,  a  value  of 
24  g's  was  obtained  for  the  average  of  the  1/10  highest  bow  acceleration 
from  the  experiment  for  a  craft  with  10  degrees  deadrise  angle  as  compared 
to  11  g's  from  the  computed  results.  These  values  were  obtained  in  a  seaway 
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with  a  significant  wave  height  to  beam  ratio  of  .666  and  a  speed  length 
ratio  of  6  which  represent  an  operating  condition  far  more  severe  than  that 
in  which  the  boat  would  be  expected  to  operate  (42  knots  in  a  State  4  sea 
for  a  50  foot  craft). 

In  summary  it  appears  that  the  computer  program  can  predict  craft 
behavior  quantitatively  most  effectively  in  moderate  operating  conditions. 

In  severe  operating  conditions,  probably  beyond  that  in  which  the  craft 
would  be  expected  to  operate,  the  computed  vertical  accelerations  which 
include  impacts  are  roughly  one  half  the  experiment  values;  however,  despite 
the  deficiency  of  the  computer  program  in  predicting  quantitative  results, 
the  predictions  are  still  indicative  of  gross  trends. 

Additional  work  is  required  to  improve  the  prediction  of  the  impact 
accelerations  especially  during  the  more  severe  conditions.  Towards  this 
end,  it  would  be  desirable  to  compare  the  time  histories  of  experimental 
motion  and  acceleration  data  with  the  corresponding  computed  values.  It 
may  be  possible  to  modify  the  hydrodynamic  coefficients  in  the  mathematical 
model  on  the  basis  of  experimental  results,  specifically  those  affecting 
impact  acceleration,  and  greatly  improve  the  correlation.  It  is 
recommended  that  if  additional  model  experiments  are  conducted,  time 
histories  of  the  motions  and  accelerations  along  with  the  wave  height 
be  made  available  for  the  above  studies. 
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WAVE  AMPLITUDE  REPRESENTING 
DISCRETE  SPECTRUM 
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Figure  4  -  Sample  of  Computed  Bow  and  CG  Acceleration 
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Figure  8  -  Comparison  of  Computed  and  Experimental 

Pitch  Variation  with  Significant  Wave  Height 


Figure  9  -  Comparison  of  Computed  and  Experimental  Heave 
Variation  with  Significant  Wave  Height 
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Figure  11  -  Comparison  oftamputed  and  Experimental  Pitch  Variation  with 


OEAORISE  ANGLE  <*•*•«•) 


Figure  12  -  Comparison  of  Computed  and  Experimental  Heave  Variation 
with  Deadrise  Angle 


Figure  13  -  Comparison  of  Computed  and  Experimental  Bow  and  CG 
Acceleration  Variation  with  Deadrise  Angle 
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PROGRAM  MA IN (INPUT. OUTPUT. TAp£5®iNP0T,TAPE6»0uTPUT.TAPE3a8l2« 

MAIN 

2 

.  TAPE2®S12»TAPe*.*5l2*TA>JE9) 

MAIN 

J 

c 

MAIN 

4 

real  it.m.mmax.n.no.nl.nar 

MAIN 

5 

INTEGER  t'NO 

MAIN 

6 

c 

MAIN 

7 

DIMENSION  y (6) .FX (2,2000) 

MAIN 

8 

c 

MAIN 

9 

COMMON  /CONST/  NCO.ECG.PI.OHR.»PU.GMAVTT.RHt),NUM.MA|120) .CD.TA, 

MAIN 

10 

.  H(  120)  .BtTA.rt*  (I/O)  .TZ.UP4G.8.XD.T.XP.M.IT. 

MAIN 

11 

1  OELT  AS.TX.ES T  <1201 . KAP.mmAX  U20T.TEST ( 120  r  > 

MAIN 

12 

.  N ( 120 | .PHALF 

MAIN 

13 

COMMON  /SNIP/  MASb.ClNT.0A.CE.CE2.Cej.0MU.E0MU.E20rtU.E30MU.BF. 6mm 

.main 

14 

.  NL.FL.IA.EU20> 

main 

15 

COMMON  /IN/  hm<120)  ,B1  U2U)  .VtLIN 

main 

16 

common/out/npmint.nplot.eno 

main 

17 

COMMON/ TERMS/ I l.T 2. T3.T4.T5.Tb.T7.T8 

main 

18 

COMMON  /TPaNS/  ST APT .RISE .PAMP 

main 

13 

common  /inter/  ii.ktt  uo>  «oife  uo> 

Main 

20 

common  /IN?/  NGU20)  ,XA«Xt.MMAX,HMlN,A(b>  .EPSElfe) 

MAIN 

21 

COMMON  /ACCEL  /  XaCCL.8aACL.CGACL.8L 

MAIN 

22 

c 

main 

23 

CALL  INPUT 

MAIN 

24 

c 

main 

25 

c 

COMPUTE  INTEGRATION  INTERVAL  INFORMATION 

main 

26 

c 

MAIN 

27 

NLESS  »  NUM-I 

MAIN 

28 

I  a  1 

MAIN 

29 

11  =  1 

MAIN 

30 

DIFFER  a  EST <  I  ♦  1 ) -EST ( I > 

main 

31 

ATT  < 1 1>  *  1 

MAIN 

32 

OIfF(II)  =  differ 

MAIN 

33 

U0  2S  1*2. NLESS 

Main 

34 

differ®  esTU.n-tSTm 

main 

35 

NTT  •  1 1  •  s  KTTUH.1 

Main 

36 

IF(0lFFER.NE«0lFF  (ID)  GO  TO  2h 

main 

37 

GO  TO  25 

MAIN 

38 

z<* 

II  ®  1I» 1 

main 

39 

KTTU1)  *  l 

main 

40 

OIfF(II)  ®  differ 

main 

41 

2» 

CONTINUE 

main 

42 

KTTUII  =  NTT  ( 1 1)  *1 

main 

43 

c 

•  • 

*  *  *  CHECN  IF  num«eR  OF  INTERVALS  EXCEEDS  DIMENSION 

main 

44 

IF  (II.GT.10)  »RI  Tt (b.28 )  (NTT  1 1 ) .D1FF ( I ) . 1  =  1 . 1 1 > 

main 

45 

IF(II.GT.IO)  STOP  4 

main 

46 

c 

•  • 

*  *  •  POINT  AT  AMICM  MULTIPLE  PUNS  START 

MAIN 

47 

0 

CONTINUE 

MAIN 

48 

call  seaaay 

main 

49 

CALL  TABLE 

MAIN 

50 

IlMEsXA 

MAIN 

51 

NOijNTal 

MAIN 

52 

eno®£no- 1 

MAIN 

53 

•  RITE (6  *  39> 

main 

54 

3* 

FORMAT  C 1H1 > 

main 

55 

c 

•  • 

••••••  Rt AO  IN  INITIAL  CONDITIONS 

main 

56 

c 

X(l>  *  VELOCITY.  X(2)  *  Z  DoT  »  X<J>  *  THET4  DOT 

main 

57 

c 

X (4)  a  X,  X ( 5 )  *  Z»  X<6>  =  THETA 

MAIN 

59 

c 

THETA  IS  -EAU  IN  DEGREES  THEN  COmVEPTEO  TO  RADIANS  IN  PROGRAM 

main 

59 

c 

MAIN 

60 

«EaO(5»10MXU>.UI,M  MAIN 

c  ^  r.  UATA  •  ust0  IN  H*M><  ►unction,  to  town  on  WAVE  main 

«eao<‘»,to)<?tant,nise  main 

10  FORMAT (8F10. 4)  "?{* 

c  ••••••••  yhl T£  out  THE  INPUT  VmLUES  main 

WHITE (6, 19)  START. RISE. KAR  MAIN 

H  FORMAT («  STAHT  *  F10.  •»,/,»  WISE  »  M»FlQ.4»/»»  KAR  a  H, FIONA  IN 

c  **4’  main 

C  THE  IS  THE  tine  AT  WHICH  TnE  INTEGRATION  INTERVAL  IS  MAIN 

C  TV)  HE  Changed 

c  HNX  IS  THE  NEW  MAXIMUM  INTERVAL  SIZE  AFTER  TIME  THE  MAIN 

C  HNN  IS  THE  NEW  MINIMUM  INTERVAL  SiZt  FOR  KUTMER  TO  SUR-OIVIOE  MAIN 

C  The  MAXIMUM  Interval  up  to 

C  IF  THIS  OR T 1 ON  IS  NOT  JStD  SET  THE  TO  THE  STOP  TIME  OF  THE  RUN  MAIN 

R£AD(S.l0)  TME , MMX , HMN  MAIN 

wRIT£<6«11)  Tme ,hmxx,HMX,hmIn,hmn  main 

/I  nHc  5,1t  N**1>,U*  INTERVAL  S12E  FOR  INTEGRATImaIN 

•ON  WILL  rifc  CMANoEU  From  «,F  10.4, •  TO  *,FI0.>»,/,  MAIN 

••  ANO  The  minimum  SIZE  FOR  HALVING  CHANGES  FROM  • ,r l o .4 .  MAIN 

•  •  TO  • ,F 10.4) 

£  ?HJU5TcTH£  !lMt  F0R  ChAN<*  UF  INTEGRATION  INTERVAL  MAIN 

C  for  CHECK  against  TIME  In  Tme  1NTEGART10N  LOOP  main 

TM  =  TME- (HMAX/2.  I  hmn 

c  MAIN 

IPT  a  0  M*,N 

IF(TME.EO.XE)  IPT  *  1  "*{* 

WEaO(S,10)  PERCnT 
XACCL  *  ECG-PEhCNT*BL 

wRlTE(6,12)  PtnLNT, XACCL  " 

**  ^R2~T<*  ™E  *  UbtD  F0Q  rHt  t,0“  MNt)  CG  ACCELERATION  COMPUTATIONS  main 
EQUAL  TO  ECO-*, F 10.4, 7H«hL  Or  ,F10.4)  MAIN 

WRITE  (6,23) 

WRITE(6,-?) 

23  FORMAT  (1H  ,//)  "’*{N 

47  FORMAT STATION  N0.m,3X,«0£a0  Rl St«,8X ,»EST»,8X,«N0«,  main 

•  10X,«fl£A<»)  mm  i  p, 

WRITE»6.SS)  ((l,HtTA,EST«I>«NO<I),bM(I)),IaJ,NuM, 

5s  F  ORMAT (8X,I2,bx,F10*M,NX,F10.4,4X,FI0.4,3x,F10.4)  M»In 

WRI  TE  (6.23) 

WRITE  (6,5*0  (  X  1 1 )  *1  =  1,6)  ”7, 

S6  F  ORMAT ( "  X  VALUtS»,4X,6(F10.4,2x))  I 

. . CHunOL  INPUT  ^  ROM  UEGrEES  TO  RADIANS  MAIN  : 

X<3>  =  X(Jt*RPo 

X<6)  -  X  (6)  *RPU  mA{JJ  j 

:••••••••  nriie  out  compuieu  arrays  ! 

•RITF(6,67)M,lT,Prl»LF,PI,0PAVlY 

1F(NPRINT.LT.4)  uO  to  62  *?{"  J 

•RITE  (6.5A)  (E ( I ) , 1=1 ,NUM)  "win  } 

WRITE  (6 ,89)  (.x(I)  ,I  =  1,NUM)  ‘ 

wRITE  (6,6L)  (PMAXI1) ,l*l.NOM)  "in  { 

WRITE  (6,661  ( Its? ( I) ,I»1,NUM)  MAIN  1 


.IS  EQUAL  TO  EC0-».F 10.4,7H«hL  Or  .F10.4) 

WRITE (6.23) 

WRI T£ (6,«7 ) 

23  FORMAT (1h  ,//) 

47  FORMAT (••  STATION  NO.«,3X.«OEaO  RiSt",8X,»EST»,8X,«N0«. 
•  iox,«a£A<») 

WRITE (6.56)  <(I*BtTA,ES?(l).NO(I),BM(I)),I=j,NyM) 

5a  F ORMAT (6X* 12. bx  , FIO, 4, MX. F10, 4. MX. F 10. 4, 3x. FlU. 4) 

WRI TE (6,23) 

WRITE  (6,5*.)  <x(l)  *1  =  1,6) 

56  FORMAT!"  X  VAL0bS"«4X ,6 (F1J,4,2X)) 

»•••••••  CHunOL  INPUT  ^  ROM  DEGREES  TO  RADIANS 

X (3)  =  X(Jt*RPo 
X (6 )  =  X (6) *RPD 


••••••  WRITE  OUT  COMPUIEU  ARRAYS 

•RITE (6,67 )M, IT  .Prt»LF,PI,G9AVlY 
IF  (NPrInT  .LT  .4)  uo  JO  62 
WRITE  (6«5<M  (E ( I ) , 1=1, NUM) 

WRITE  (6,8-))  (.,  ( I )  ,  I  =  J  ,NOM) 

wRITE  (6,64)  (MMAX(l)  ,l*l»MOM) 

WRITE  (6,661  ( IfcSl ( I) ,twl,NUM) 


N  61 
N  62 
N  63 
N  64 
N  65 
N  6b 
N  67 
N  68 
N  6V 
N  70 
N  71 
N  72 
N  73 
M  74 

75 

76 

77 

78 

79 

80 
81 
82 

83 

84 

85 

86 

87 

88 

89 

90 

91 

92 

93 

94 

95 

96 

97 

98 

99 
100 
101 
102 

103 

104 

105 

106 

107 

108 

109 

110 
111 
112 

113 

114 
lib 
116 

117 

118 
119 
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62  CONTINUE 

WRITE  16.28) (KT7 (1) .DIFF1I) .1=1.11) 

28  FORMAT!*  KTT.U1FF  *, I 10.2X ,F10.4) 

S I  FORMAT  <4H  M3  ,F 10.4.4H  1*  .F10.4.11H  Pl*RH0/2*  .F10.4, 

1  SH  PI»  .F 10.4.1 ON  GRAOlTTw  .F10.4) 

58  FORMAT  F.  <1)»,10F  10.4) 

5 9  FORMAT  N  ( 1 )  " .  1  OF  1 0  .4 ) 

64  FORMAT  1“  MMA*U)»,10F10.4) 

66  format  i»  test (1)".iofio.4> 

18  *  1 

1 PR I NT  *  NRrImT 
WRITE  14.91 ) 

<;••••••••  WRITE  HEADINGS  AND  CONDITIONS  AT  TINE  •  0. 

91  FORMAT  (lHl.2X.,'nMtM,9X,"X0OT”, 9*. •■2D0T“.9X. "THETA  D0T«.6X. 
.  lMX.9x,lHZ«9x,i»rllH£TA*4X,2HNC«9A.2HFL» 

.  4X.8H80W  ACCL.4X.7HCG  ACCL,//) 

WRITE (4. 92)  Tint, (X ( 1) .1*1.6) .NL.FL.bw ACL. CgACL 
WRITE (9)  TIME. (X( 1) , 1=4,6) .HWACL.CGACL 
NOlfNT  »  KOuNT.l 
FX<l.l8)=X<5> 

Fx 12. 18) *X (6) 

1KUTM  =  (Tm£-XA) /hmax  ♦  <XE«-TME)/HMX  ♦  .05 
F1RST«0.0 
NEf)S*b 
IKUTS*0 

START  OF  INTEGRATION  loop 
8S1  continue 

NPrINT  *  IRRIM 
C***«***«  CHtCK  PITCH  ,GT.  .5236  RADIANS 
IF (X(6) ,GT.  .5«i J6)  GO  TO  853 
C********  PERFORM  INTEGRATIONS 

IF, TIME. LT.TM. UR. THE. EO.Xt)  GO  TO  98 
IF(IPT.EO.I)  GO  TO  98 
HMIN  *  HMN 
HMAX  =  HMx 

first  =  0.0 

98  CONTINUE 

CALL  KUTMEP  (Nt 'jS.  T  IME.  HMAX.  X.tPSt.  A,  HMIN, FIRST) 

1KUTS=I*UTS«1 

IF (F IRST . t9, 2) GO  TO  861 

IF (KOUNT .nE. 1 . anD.AOUNT .Nt.4l )  GO  TO  99 
wRj  TE 14,91 ) 

IXOUNT=l 

C  ••••*•*  4  wnlTt  OUT  TIME  INTERVAL  RESULTS 

99  WRITE (4,92)  TIME. (X(I) ,1*1,6) .NL.FL.bwACL.CliACL 
WR1TE16,9J)T1»I2.T3,T4,T5,T6,T7,T8,8MM,8F 
WRITE (9)  TIME. (XI 1 )» 1*4*6) .awACL.CuACL 
IF(TIME.LT.TM.uR.TME.EOtXt)  Go  TO  200 
IF(IPT.EO.l)  (jO  TO  20u 

call  ploteq  (Fx.xa.hmax.ig.ipT) 

1PT  *  1 
18  *  0 
XA  »  TIME 
F lRST  *  0.0 
HMIN  *  HMN 
HMAX  =  HHX 
200  CONTINUE 
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MAIN  179 

FX(l.lB)=X(S> 

*r  MAIN  180 

Fx<2, lb>=A<6> 

^ MAIN  181 

93 

FORMAT I "  ", 10E1U.4) 

▼  main  182 

92 

FORMAT  I IX. 11 (Flo. 4. 2X)»  . 

main  193 

100 

CONTINUE  r 

MAIN  194 

nount*kount*i 

MAIN  185 

21 

CONTINUE 

MAIN  186 

IFiTIHE.LE.XE.anO. 1KUTS.LT. 1KUTm)O0  TO  651 

MAIN  187 

WRITE  (2. 63/*) 

MAIN  18B 

85<* 

CONTINUE 

MAIN  199 

AS2 

FORMAT!  ••  END  OF  KUTMER") 

MAIN  190 

853 

CONTINUE 

MAIN  191 

CALL  PLOTEo  (FX.XA.HMAXfIb.lPT) 

MAIN  192 

c 

•  • 

•  *  *  •  •  CHtCiv  FOR  LAST  run  if  not  CYCLE  -<ack  to  reao 

MAIN  193 

c 

new  DATA  FOP  we at  RUN 

MAIN  194 

IF(En0.nE.1)00  to  8 

MAIN  195 

00  TO  994 

MAIN  196 

c 

•  • 

*  •  KUTMEP  ERPOP  MESSAGES 

MAIN  197 

961 

WRITE  (6»bo2> 

MAIN  199 

862 

format («  error  criterion  in  kutmer  can  not  be  met**) 

MAIN  19v 

WRITE  (6.56)  (XII). 1=1, 6) 

MAIN  200 

■RITE  (6.6M  TIME 

MAIN  201 

66 

FORMAT  ("  TIMt  F  10.4) 

MAIN  202 

IF (EnO.nE. 1 ) GO  TO  8 

MAIN  203 

uO  TO  85J 

MAIN  204 

99  V 

continue 

MAIN  205 

ENO  FILE  9 

MAIN  206 

ENO 

MAIN  207 

SUBROUTINE  PLO 1 2 (F .fmin.fmax.nvar.nfun.ni.n.xo.uelx) 

PLOT 2  2 

c 

PL0T2  3 

c 

PLOT  FIRST  N  POINT 5  OF  UP  TO  26  FUNCTIONS  FIX) 

PL0T2  4 

c 

F  ( I , J)  CONTAINS  Trt  VALUE  FOR  1HE  jTr*  POINT  OF  THE  ITM  FUNCTION 

PLOT 2  5 

c 

FHIN(I)  ANO  HAX(l)  CONTAIN  THE  M£N  ANO  MAX  ORDINATE  VALUES  FOR 

PLOT2  6 

c 

the  itm  function. 

PLOT2  7 

c 

NVAR(I)  AN  ARRAY  OF  TITLES  FOR  TmE  VARIOUS  FUNCTIONS 

PL0T2  8 

c 

Tu  Bt  PLOTTEO  AGAINST  THE  ABSCISSA 

PLOT2  9 

c 

NFUN  NUMBER  OF  FUNCTIONS  TO  BE  PLOTTED  -  DIMENSION  OF 

PL0T2  10 

c 

NVAR,  FMIN.  FMAX 

PLOT2  11 

c 

NI  UStO  ONLY  IN  F(Nl.l)  AS  PASSEO  DIMENSION 

PL0T2  12 

c 

N  NUMBER  OF  POINTS  IN  «  SINGLE  PLOT  FPAME 

PL0T2  13 

c 

XO  FIrST  ABSCISSA  VALUE 

PL0T2  1 4 

c 

oelx  mbscissa  increment 

PLOT2  15 

c 

PLOT2  16 

U I mens I ON  JSTEP (2o) ,F(N1,N) ,FMIN(NFUN) .FMAX (NFUN) »VLAST(26) . 

PL0T2  17 

1 

Vf I »ST (C6) ,HEA0{6> .STEP (26) 

PLOT2  18 

Integer  c-h26)  »nvar<  nfun> ,uot,«ster,plos»3Lan* 

PLOT2  19 

integer  c 

PL0T2  20 

INTEGER  4(101) 

PLOT2  2) 

c 

PL0T2  22 

data  bL ANi « OU T  » AS 1ER .PLUS/ 1h  , lH., IH*. lH./ 

PLOT 2  23 

DATA  Ch(l) ,CM(F) ,CH(3) ,CH(4| ,CH(5) ,CM(6) ,CM(T) ,CM(8) ,CH(9| ,CH(10)  PL0T2  24 

l 

>  /  IMA  .  1  FIR  ,  1HC  ,  HO  ,  1  HE  ,  1  hF  ,  1  1&  ,  1  FtH  .  1HI  , I HJ 

/  uLOT2  25 

OAT A  CH( 11) ,Cn( 12) , CH (13) .Cm (14) .CM (15) .CM (16) ,CM< 17) ,CH< 16) 

PL0T2  26 

l 

|  /  IrtK  ,  IF1L  .  1HM  ,  1HN  ,  IHO  .  IMP  .  iHO  .  IHR/ 

PLOT2  27 

DATA  CHI  14) .cm (20) .CH(21) «C*(22) ,CM(2J) ,CH(24) .CM (25) .CM (26) 

PLOT2  28 

2 

!  /  IMS  ,  1  NT  ,  1MU  •  1 MV  ,  1 MW  ,  ltx  .  1  MY  ,  1HZ  / 

PL0T2  29 

c 

PLO I 2  30 

IF (NFUN.LE.O.UP.N.LE.O)  RETURN 

"L0T2  31 

36 
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C  PRINT  HEADINGS. 

4RI  T£  (6««6) 

46  FORMAT  (///I 
00  40  I=1.NFUN 

30  TEMM*48S<FMAX<1)-FmIN<I)) 

ExP=l. 


If  ITENH.EQ.W.)  GO  TO  2 
C  BRING  TENM  TO  A  VALUE  BETWEEN  1  ANO  10 
IF < TENH.LT. 1 • I  GO  TO  1 
3  IFTTENM.LT.10.)  GO  TO  2 

CxP=EXP«l(i. 

T£NMsT£NM».l 
GO  TO  3 

1  £XP»£XP*.l 


T£NMaT£<xM»10. 

IFTTENM.GT.l.)  GO  TO  2 
GO  TO  1 

C  SET  UP  VALUE  HETutEN  GRIO  LINES.  RSTtP. 

2  PST£P*5. 

lF<TENM.G£.b.)PST£P*10. 

IFTTENM.LT.2.)PSTtP=2. 

5  RSTEP  < I ) =uSTEp*E  XP* . 1 

C  COMPUTE  valul  of  starting  line,  vfiwst. 
FIRST»FMlN<I|/xStfcP<I) 

IFTFMlNTl) ,LT.0.)F I«ST*FIPST-1. 
FIRST*AInt (F1«ST) 

VF I RST ( I ) aF I wa  T*PST£P ( 1 1 
C  CMECK  ENO  LINE  VALUE. VLAST. 

VLAST(I)=VF1PSI <i>  .10.*RSTtPJI> 

IF  T  VLAST ( I ) .01 .Fha X  <  I ) ) GO  TO  4 
C  IF  GRAPH  IS  TOO  SMALL  TAKE  NEXT  LARGER  STEP. 

aa^pstep 

IFTAA.LT. S.|PSTEP=5. 

IF ( A A .£0.5. ) PSTEP* 10. 

IFTAA.LT. 10.1  GO  TO  5 
PSTEP«2. 

EXP310.*t<P 
60  TO  S 

C  COMPUTE  VALUE  rt£T*tEN  POINTS. STEP. 

4  STEP»I)sPSTtPU>#.l 


Rp»0. 

UO  6  kkaI.a 

rlEAO  OXK)  *VF  ILST  T  l )  *2,*RK«RSTtP  T I  > 

6  »K«RK*l. 

40  PRlTt  (o.-»S)  CHTI),  NVARTI),  THE  AO  TKM  »KK=  1.6) 
45  FORMAT ( 1 X  »  a  1 . jn  »  . A1 0.5X. 1PE12.4 ,5 (8X. 1PE1 2.4) ) 
O0  3U  j: I. 101 


* 


A  ( J)  *BLAN8 

IE  (MOL)  (J.lO).Eu.l)  ATJ)«00T 
50  CONTINUE 

«RlT£<fr,55)  A. A 

55  Format  <2SX.lvlAl/lSX»4HTlM£*6x.l0lAl) 
C  PLOT  Each  point 
UO  100  J=1.N 


b*xO .FLO A  T ( J- 1 ) »OtLX 
UO  1 0  X« 1.101 
ATK) “PLAN* 

IF  (MOO(K.l'J)  . E 0 . 1 )  A(K)>0OT 
IFJMOOTJO  l.t'i.I)  A(K)*DJT 


PL0T2  32 
PLOT2  33 
PL0T2  34 
PL0T2  35 
PL0T2  36 
PL0T2  37 
PLOT 2  38 
PLOT2  3R 
PL0T2  40 
PL0T2  41 
PLOT2  42 
PLOT2  43 
PL0T2  44 
PL0T2  45 
PLOT?  46 
PL0T2  47 
PL0T2  48 
PLOT2  4R 
PLOT 2  50 
PLOT2  51 
PL0T2  52 
PLOT?  53 
PL0T2  54 
PLOT?  55 
PLOT?  56 
PLOT?  57 
PLOT?  58 
PLOT2  SR 
PLOT2  60 
PLOT2  61 
PLOT?  62 
PLOT?  63 
PLOT?  64 
PL0T2  65 
PLOT 2  66 
PLOT2  67 
PLOT?  68 
PLOT?  6G 
PLOT2  TO 
PLOT?  71 
PLOT?  72 
PLOT2  73 
PLOT?  74 
PLOT?  75 
PLOT?  76 
PL0T2  77 
PLOT 2  78 
PLOT?  7R 
PLOT2  80 
PLOT 2  81 
PL0T2  82 
PLOT?  83 
PL0T2  84 
PL0T2  85 
PLOT2  86 
PLOT?  87 
PLOT?  88 
PLOT?  8R 
PLOT2  RO 


37 


onoononnnnoonnon 


TO  CONTINUE 

00  00  lsl«NFU'« 

L0Ca( (F<I,j)-vF1RST(I>  >/STEP(l  )*1.S) 

CsA(LOC) 

A (LOC) =CH ( I ) 

IF (C.NE. BC*NK. ANO.C.NE. DOT)  A(L0C> “ASTER 
80  CONTI Nut 

IF(N00I  J.10)  .EU.1IGO  TO  95 
BRlTE<6.85>  A 
85  FORMAT  (25X* 1 U1A1 ) 

00  TO  100 

99  aRtr£(6tlb)9«4 
1 o  FORMAT  (12X*1hE1Z.4.IX|101A1) 

100  CONTINUE 
RE TORN 
CNO 

SUBROUTINE  MUTMER(NO,T«H»VO.EPSE.A«HCX»FIRST) 

DIMENSION  vO<6>  «Y1 (6) «Y2(b) .F0(6) «F1 (6) .F2(S) ,EPSE(6i .A (61 

COmMOn/OuT /NPR INT  .NPLOT  «£N0 

COMMON  /ACCEL  /  xaCCL.BmACL.C&ACL.BL 

DATA  NAMI.nAME  /2HVU2HY2  / 


PL0T2  91 
PL0T2  92 
PL0T2  93 
FL0T2  94 
PL0T2  95 
PL0T2  96 
PL0T2  9 7 
PL0T2  9ft 
PLOT  2  99 
RL0T2100 
RL0T2101 
PLOT2102 
HLOT2103 
RLOT2104 
PLOT210S 
PLOT2106 

MUTMEW  2 
AUTMEH  3 
KUTMER  4 
KUTMER  5 
KUTMER  6 
MUTME*  T 


NO  a  NUMBER  OF  EQUATIONS.  NO.  OF  COMPONENTS  OF  VO 
T  *  INDEPENDENT  VARIABLE 

H  a  INCREMENT  FOR  BHlCM  SOLUTION  IS  TO  BE  RETURNED  ♦  OR  - 
VO  «  THE  VECTOR  UF  CJEPtMuENT  VARIABLES.  ENTER  BllM  INITIAL 
VALDES  AT  T  AND  RtTuRN  aITh  VALUES  AT  T*M 
EPSE  a  RELATIVE  ERROR  CRITERION  FOR  COMPONENTS  of  YO  .GT  ABS(A) 

A  a  ABSOLUTE  ERROR  SR  1 TERlON  ►OR  COMPONENTS  OF  VO  .LT.  A0S(A> 
NOTt—  EPSE  ANO  a  must  BE  SPtClFltO  FOR  EACH  COMPONENT  OF  THE  SYSTEM 
HCX  a  THF  SMALLEST  STEP  Sl/t  USEU  IN  THE  INTEGRATION 
first  SHOULD  HE  o  «H£N  KuTMER  IS  cNTERtO  FOP  the  first  time 
after  That  FIRSI  IS  1  IF  MOTHER  IS  ENTERED  MITH  THE  SAME  h  OR 
IF  IT  IS  ENTtREU  WITH  A  CmauGEO  H 
IF  FIRST  IS  2  THE  OR0R  SRITEPIA  CANNOT  bE  MEET  AND  The  STEP  SIZE 
REOUCEO  TO  R/12B. 

IF 


10 


(FIRST)  20.lu.E0 

-  -  -  -  -  FIRST  ENTRY 

nC  =  H 


1PLOC 

FIpST 


1 

l. 


...... - -  OTHER  ENTRY 

20  LOC  =  0 
HCX  a  HC 

IF  (HC.NE.O.I  GO  10  30 
«RITE'6.dU>)) 

800  FORMAT (Sa.hSHMLTmER  ENTERED  WITH  ZERO  INTEGRATION  INTERVAL  ) 
FIRST  =  C.. 

RETURN 

..........  b  CALLS  TO  UAU* 

30  CALL  OAU* (T.YO.FO) 

IF (NPR INT. EO.^I "MITE <6. 400) YO.T.FO 
400  FORMAT  <6<C<.F lu.4) . 4HI I ME  *2A.F10.4| 

IF (NPRINT.EO.S) *R11E (6.400) HC 

39  00  40  I a 1 • .0 

40  Y 1 ( I )  *  Y 0 ( I > ♦ ( HC /).)»FOtl) 

IF(NPRlNT.r  Gi.si  WRITE  (6.400)  Yl»T 


KUTMER  8 
KUTMER  9 
KUTMERIO 
MUTMEMII 
MUTMEH12 
MOTHER 13 
MUTMER14 
MUTMER1S 
MUTMER16 
MUTMER17 
MOTHER  18 
MIJTMER19 
IMJTMER20 
MUTMER21 
MUTMER22 
MUTMER23 
MUTMER24 
MIJTMER25 
MUTMEH26 
MIITMER27 
MUTMER28 
MUTMER29 
MUTMER3P 
MOTHER 3 l 
MOTMpR 32 
MUTMER33 
MUTMER74 
MUTME»35 
MUTPER36 
MUTME Ri7 
MUTMER38 
MUTM£k39 
MUTME.40 
MUTME -<41 
MUTME R42 
MUTMER4J 
MUTMER44 
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CALL  OAUX  <  T .HC/3, , V 1  ,Fl  > 

If (NBRINT.* 0.3) «R1T£ <6. 400) FI.  T 
UO  50  1*1, ,0 

V 1(1)  «  TO ( 1 ) ♦ (HC/8. ) *F0 ( 1 ) * (HC/6. ) *F1 ( I > 
IF(NPHINT.F.0.5)»R1TE(6,40U>  tl.T 

CALL  OAUX<T.HC/J..Vl,f 1) 
lf(NPRlNT.EQ.3>«RnE(6.400>F  l.T 
00  60  Ial.NO 

tid)  »  to<i>*(iic/a.)*fO(i>».j/5«Hc*fi(n 
if  (NPRlNT.FQ.S)  »RIT£  (6,400  >  Tl.V’ 

CALL-  OAUX  (T.HC/2.»Y1,F2) 

If (NPRlNT.EO.b) 00 11E (6.400  IF  it T 
00  70  Ial.NO 

till)  *  rU(I).<HC/2.)#F0(I)-l.5*HC*fl(I>»2.*NC»F2(I> 
IflNPRlNT.EQ.b) aR1T£ (6.400) VI, T 

CALC  OAUX(T.HC.Vl.F 1) 

If  (NRRlNT.FO.aUHllE  (6.400)7  l*  I 
UO  00  I«1.M0 

72(1)  ■  VO ( I ) «hC/0* *f 0 (I) ♦ <2./3. ) •HC*f 2 C I ) ♦ < HC/6. ) *f 1 ( I ) 
If  (NPRlNT.FO.  bllUME  (6*400)  Y2*T 
iNr  »  0 

-  -  -  - - CHECK  ERROR  CRITtRU 


«  w 

UO  110  I 

a  1  .NO 

ZZZ  »  ABS(tlU) 

IF  (ZZZ) 

85 

Error  * 

AdS ( .2* 

If  (ERROr-A(D) 

87 

tRROR  a 

ABS (.2- 

IF (ERROR 

•  m  •  • 

-tRStll 

m  • 

90 

X  a  128. 

*A6S(HC 

IF  (X)  91 

,95,9b 

91 

■RrTE (6, 

92)  I.T.i 

ABSOLUTE  ERROR 


RELATIVE  tRROR 


SINCE  ERROR  .l»T .  ERROR  CRITERIA  CHECK  If  HC.GT 
IF  VES  THEN  HALVt  INTERVAL.  OTHERWISE  STOP. 


ERROR  TOO  LARGE 


i  t ORMAT  < /  IBi-t  FOR  EOUATION  NO.  I2,2I«.  TM£  RELATIVE  ERROR  AT  T  * 
.  £15.6.  4h  IS  ,E 15 .8, 1 JH  STEP  SIZE  *  ,Elb.8> 

FIRST  a  2. 

RETURN 

.........  HALVE  INTERVAL 

»  HC  *  HC/2. 

IPLOC  a  2* 1PLOC 
i-OC  »  2*LUC 
HCX  *  HC 

6«ITE(2,710)T,I,ERROR,HC 

)  FO»RAT(/6H  TI«t  «  *F10.3.5X»26HHALVE  INTERVAL.  EQUATION  .13, 
4 13m  HAS  t»»OR  =  .E16.8.6X.1 7H  STtP  SIZE  NOR  «  ,£15.8) 

6RITE  <2.  72o)  namZ,  (Y2(  J)  ,J=1  .inO) 

»RI TE  <2, 720)  NA«1, (tl<J),J»l,NO) 

)  FOrmaT(  2X.A2  /  3 < 10E13.5/1 ) 

00  TO  30 

-  TtST  If  INTERVAL  LENGTH  CAN  BE  DOUBLED 

)  IF  (ERR0R**>4.-tRSt  ( I ) )  llU.ilU.101 
l  INC  =  1 


KUTMER4S 

KUTMER46 

KUTHER47 

KUTMER48 

KUTMER4V 

KUTMERSO 

KUTNER51 

KUTHER52 

KUTMER53 

KUTHER54 

KUTNER55 

KUTHER56 

KUTHER57 

KUTHER58 

KUTMERS9 

KUTMER60 

KUTHER61 

KUTME°62 

KUTHER63 

KUTHER64 

KUTNEP6S 

KUTMER66 

KUTHER67 

KUTHE»68 

KUTMER69 

KUTMER70 

KUTRER71 

KUTRER72 

KUTHER73 

KUTHES74 
KUTHER75 
KUTMER76 
KUTMfR77 
KUTHER78 
.H/KUTHEH79 
KUTHER80 
KUTHER81 
KUTMER82 
KUTMER93 
KUTMER94 
,  KUTMER95 
KUTMER8A 
KUTMERH7 
KUTMER88 
*UTM£Ra9 
KUTHER90 
KUTMER91 
KUTMER92 
KUTHER93 
KUTHER94 
KUTMER95 
KUTH£R9t> 
KUTMER97 
KUTMER98 
KUTMEW99 
KUTHEldO 
KUTME101 
KUTM£ 1 02 
KUTHEI03 


o  u  uuu 


110  CONTINUE  KUTHE104 

c  -----  -  —  -  -  update  t  ano  solution  autmeios 

111  T  *  T *HC  KUTME106 

00  112  1*1, NO  KUTHE107 

112  VO(I>  ■  V2<n  KUTHE108 

C  -----  -  -  otT  SOLUTION  IN  NtXT  INTERVAL  *UTME109 

LOc  *  LOO  l  KUTHE110 

IF  ILOC-IHLOC*  120*210*210  *UTMElll 

120  1F(1NC)210.130*21J  RUTHE112 

130  IF  (LUC-(L0C/2)»2)  210*140*210  *UTME113 

140  IF  ( 1RLOC-1 ) 21 0 *21 v* 200  KUTME114 

c------ - -  uouble  interval  length  kutmeiis 

200  HC  *  2.*rtC  RUTME116 

LOC  *  LOC  /2  RUTHE11T 

IPLOC  *  IHLOC/2  RUTHE118 

210  IF ( IPLOC-LOC)  30,320,30  KUTME119 

320  BwaCL  *  FO  (2) -XACCL»F0(3)  RUTME120 

CGACL  *  F 0 ( 21  RUTHE121 

RETURN  KUTMF122 

ENO  RUTME123 

SUBROUTINE  OAu* <T1ME*X*RHS>  OAUX  2 

UAUX  3 

TIME  TIME  AT  rtHICH  SYSTEM  IS  TO  BE  EVALUATED  UAUX  4 

X  STATE  VECTOR  UAUX  5 

RHS  the  RIGHT  HAND  SiUE  OF  THE  EQUATION  S  *  F  A  OAUX  6 

UAUX  7 

REAL  KAR  OAUX  8 

REAL  IA, 1T*M,K.ma,maSS*NL«N,MMax  OAUX  V 

U I MENS I ON  X 1 6) « RHS (6) »F (3*1) *4(3*3) • INDEX ( 3*3) *  UAUX  10 

1  -MlcO) ,U<120) ,V<120) ,wOT1ME<20)  OAUX  11 

C  OAilX  12 

COMMON  /SHIR/  MASS. CINT,QA,CE«CE2«CE3*0MU»EDMU»E2UHU.E3DMU*BF*BMM, UAUX  13 

,  NL.FL*IA»E(l20l  OAUX  14 

COMMON  /CONST/  NCG,ECG,RI.ORW.MPD»&HAVTr*RHU*NUM.MA»120).CD,TA,  UAUX  15 
.  rf ( 120 ) *HET  A , HV (120) *  1 Z . URAG*  N«  XO*  T • XR»M* IT,  UAUX  16 

I  JELTAS.T  x*EST  (120) »*AM,mmAX (120)*TEST|120>*  OAUX  17 

.  N ( l 20 ) , RmaLF  OAUX  18 

COMMON  /IN/  HH(120) ,B1 (1201 «VtLlN  OAUX  19 

COMMON/OUT/NRR INT  *NRLOT  »ENO  UAUX  20 

COMMON  /TRANS/  ST4RT,RIS£*RAMH  OAUX  21 

COMMON  /WAVEI/  2MA.ZMMA,EMAS,2ZtfMA*ZriEHA*Z2«MA*E2MAZ,  OAUX  22 

1  ZvOOT (120) . AS1NPT ( 120,20) , ACOSRT (120,20) *CX6*SX6  OAUX  23 

COMMON  /WAVE2/  ■0(<0)*K(20|*C(201 »RVO (20) *RVQ2(20) *RK (20) »  UAUX  24 

.  R0I20) ,RVK (20| *H*RHS(20l  UAUX  25 

C  OAUX  26 

COMMON/SINE/ROINT ( 10001  OAUX  27 

RAMR  *  RMR(TIMt*START,RlSE)  OAUX  28 

RlH  *  Rl/2.  UAUX  29 

00  S  JJ  =1.10  UAUX  30 

wOTlME(UJ)  =  «0 < JJ) •TIM£*RMS{ JJ>  OAUX  31 

KNTRAO  =  31 8,3098861  OAUX  32 

5  CONTINUE  Daux  33 

CX6  *  COS  I  * (6) )  UAUX  34 

SX6  =  SlN(X(6))  DAUX  35 

C«*«****SET  VALUES  OF  MA  AND  B  UAUX  36 

UO  75  1  =  1* NUM  OAUX  37 

RI  *  0.0  UAUX  38 

RTO  =  X(h> »E(I)»CXb*N(I)«SX6  OAUX  39 

VI  =  0.0  OAUX  40 


40 


oort  n  n  non 


r) 


00  20  J=1,10 

MT  •  PTO*K(JI*wOTlME(J) 

ANOX  »  PT*PNTMA0 

AS  *  AMOu(ANOX,20UO.)-499. 

1F(AS.LT.-499.I  AS  *  AS*200O. 

IS  *  AS 
1C  *  IS*S0O 
OX  »  AS-iS 

IF ( IS.LE.O)  IS  a  2-IS 

IF (IS. OT. 1000)  IS  a  2002-lS 

IF(IC.LE.O)  Ic  *  2-IC 

IF ( IC.6T . 1000)  IC  3  2002-IC 

OX  «  OX/PNTRAO 

XSIN  a  POINTUS) 

YCOS  a  POINT <1C) 

ASINHTd.J)  a  XSIN*OX*YC0S 
ACOSPT(I.J)  a  VC0S-0X*XSIN 
HR  a  R0( JI*ACOSPT( I.J1 
VV  a  -WO ( J) ( j) *AS INPT ( I  * J) 

«I  *  Kl»WW 
VI  a  Vl,VV 
20  CONTINUE 

H<!>  a  RI«RAMP 

>•••••••  COMPUTE  HW  SUBMERGENCE  OF  A  POINT  ANO  R  THE  HAVE 

HHUl  IS  IN  THE  FIXED  COORDINATE  SYSTEM 
HH(I»  »  X (5) -E ( 1 1 *$X6*N ( I 1 aCXP-R ( 1 1 
lF(HHCll.GT.O)  GU  TO  6S 

Cw AFT  IS  NOT  SMBMEPGtO 
MA(I)  *  0. 

01(1) *0 ■ 

0(11  a  0. 

GO  TO  75 

63  V  (  X )  a  Vl*  >AHP 

0(1)  a  MM ( I ) / (CX6-V ( I ) *$Xt>) 

0<I)  IS  IN  THE  800y  Axis  SYSTEM  ANO  IS  THE  SUBMERGENCE 
1F(0<I)  TEST  (II )  GO  TO  #0 

CP AFT  IS  PARTLY  SUBMtRGEO 
«<I>  a  0<I)*U./TA)«PIH 
61(1)  *  0(I)*(l./f4)*PlM 
MA  ( I )  =  KAR*Pm«lf*u  ( I ) *6 ( 1 1 
GO  TO  75 

C*INt  IS  IMMERSED 

61  ARRAY  IS  DStO  FOR  THE  INTEGRALS  OVER  THE  PORTION 
„  .  THE  MULL  FOR  WHICH  THE  CHINE  IS  NOT  IMMERSED 


70  MA(I) aMMAX(I) 

6< I ) *BM ( I ) 
at  (I)=o. 

75  CONTINUE 

IF(NPRInT.LT.h)  00  TO  8S 
hRITE(6«7m)TIME 

74  FORMAT  (••  TlMt  a  •*,FlU.4| 
WRITE (6*  76)  (X(I), l  =  l.o) 

WRJ  TE (6, 7  7)  ( R ( I ) *  I x 1  «NUP 1 

WRITE  (6*M>  lh«tn  .1  =  1«NJM) 
•RITE  (6.7RJ  <  6(1) *1*1 f NUM) 
•RITE  (6.80)  <Y(Il »I  =  l,N(jM> 
•RITE  (N, 61)  (C ( I )  ,  1  =  1 , N(jM) 
•RITE (6*62) (M«(l) ,lal,NUM) 

76  FORMAT!"  x  « I  >  "«6(2X,C12.61  ) 


A.  ; 


.  >  <<>  '< 

.. 

<A'.\VVk 

>VV  '  J  ' 


v  \T--- % 


OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

UAUX 

OAUX 

UAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

DAUX 

OAUX 

OAUX 

OAUX 

OAUX 

UAUX 

OAUX 

OAUX 

UAUX 

OAUX 

OAUX 

OAUX 

UAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

UAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

OAUX 

UAUX 

OAUX 

OAUX 

OAUX 

OAUX  1 

OAUX 

OAUX  1 

OAUX  1 

OAUX  1 

OAUX  1 


41 


77  F OrMAT  («  R(I)".10F  10.4) 

7*  FORMAT  (»  "tW(ll“.lOf-10.4) 

7V  FORMAT  («  H(  U"*1UF  10.4) 

8#  FORMAT  <<•  V  ( 1 1  "  *  1  OF  10  *4) 

at  format  (••  o< u»«ioF io.4i 

8 i  FORMAT(«  M4(I)  ".10F10.4) 

85  CONTINUE 

c 

{♦*»*•»»•  COmRU Tts  NL  *40  FL  AND  THE  ASSOCIATED  INTERGALS 
CALL  FUNCT(X) 

C 

IFlWRINT.LT. 4) GO  TO  17 
WRITE (6, 15)  T  X,Fl,URAG*T  Z*W»NL,xO*  T  *XR 
IS  f  0RM4T  ( ••  ",  lutl2.e>) 

1 T  CONTINUE 

C***»*«*  COMPOTE  THE  F  VtCTOn 
F  ( 1 , 1 )  s  Tx*FL*SX6-DRAG«CA6 
F 1 1 , 1 ) *0, V 

F<2»1>  *  T7*FL*Cx6»DRAG*$X6»W 
F  <  3  ,  1  >  *NL-'  )R  AG*  XD  *  T  *  XR 

if (nrkint .lt.3> go  to  ia 

•RITE (6. lui <F (1. II *1*1,31 
10  CONTINUE 

c  **«,*•  *  COMMUTE  The  a  matrix 
4(1,1)  a  M»MA5S*SX6*SX6 
*(J,2>  *  H4SS«SX6*CX6 
A<1, 3)  *  -i)A*S*t> 

A<I.2l  *  0. 

A£|.3)  *  0. 

4(2,11*4(1.2) 

A(2,2I  *  M,HASS«CXG*CX6 
A(2,3)  *  -)A*Cxt» 

A (3, i ) *A ( 1 «  3) 
a ( 3»2) *A (2.3) 

A(3»3)*1T*IA 

IF (NRRInT .LT .3)00  TU  2S 

WRITE (6, 12) <A< 1,1) *1*1,3) 

WRITE (6,13) (At  1,2) *1*1,3) 

WRITE (6,14) (A (1,3) ,1*1,3) 

C*******  INWtRT  THE  A  MATRIX 

2s  CALL  MAT1NS(A,J,3.F ,l,l,OETtRM, ID* INDEX) 

IF (I0.E0.2) WRITE (o*26> 

26  FORMAT!"  MAlRlA  IS  SINGULAR  «> 


C**«**«A  ON  RETUON  » ILL  CONTAIN  Th£  INVERSE  MATRIX 
C  10*2  MATRIX  IS  SINGULAR 

c  *i  inverse  was  fuunc 

C 

COM"UTt  THE  RIGHT  HANO  SIDE 


RHS(ll  =  F (1,1) 

RHs<2)  *  F(2,i) 

RMSTJ»>  =  F  (3*1) 

rm$ < t )  *  o.o 

Rh<;  (4 )  *  X  ( 1 ) 

RHs(SI  *  X(2) 

RHS(6)  *  X (3) 

10  FORMAT  I"  F 11*1 )  ",3(2X,E12.") ) 
12  FORMAT!"  A ( 1 , 1 )  ",3(2x*tl2.4)) 
IJ  FORMAT!"  A(U2)  »  *3( 2X  . 1 1 2 .<*> ) 
!♦  FORMAT!"  A ( 1 1 3)  ",J(2X,tl2.*0 ) 


OAUX  100 
OAUX  101 
OAUX  102 
OAUX  103 
OAUX  104 
DAUX  10S 
OAUX  106 
OAUX  107 

oaux  ioa 

OAUX  10V 
OAUX  110 
OAUX  111 
OAUX  112 
Oaux  113 
OAl)X  114 
OAUX  115 
OAUX  116 
OAUX  UT 
DAUX  118 
OAUX  11 V 
OAUX  120 
OA.JX  121 
OAUX  122 
OAUX  123 
UAUX  124 
OAUX  125 
OAUX  126 
OAUX  127 
Oaux  128 
oaux  129 
UAUX  no 
OAUX  131 
OAUX  132 
OAUX  133 
OAUX  134 
OAUX  135 
Oaux  136 
DAUX  137 
OAUX  138 
OAUX  139 
UAUX  140 
Oaux  141 
oaux  K2 

OAUX  143 
DAUX  144 
OAUX  145 
UAiJX  146 
DAUX  147 
OAUX  148 
OAUX  149 
UAUX  150 
UAUX  151 
DAUX  152 
OAUX  153 
OAUX  15" 
OAUX  155 
OAUX  156 
OAUX  157 
OAUX  158 


c  •  • 


30  IF  (NPRINT .lT .2)  Go  TO  40 
WRITE(6.12)  < A ( 1  •  1  >  ,1*1,3) 

WRITE(6,1J>  ( A  <  i  » I  « I  a  l  ,3 ) 
aRI TE (6 • 14)  (A(1,J) ,I»1»3) 

WRITE (6.3b)  (RHS(1> ,1*1.6) 

33  FORMAT  (••  RHS(l)  ",6(2X,t.l2.6)  ) 

40  CONTINUE 

RETURN 

end 

SUBROUTINE  FUmC T  <  X ) 
m£al  kah 

REAL  lA.lAA.lHAWT, API, MA, MASS. N 

integer  tN;i 

DIMENSION  IPAmT (120) ,C1  (120) .C2(120) , 

*  Ul (120) .02(120) .03(120) .04(120) .06(120) .06(120) » 

.  OPART (120) .21(120) .22(120) .23(120) .24(120) .25(120). 

.  26(120) ,27(120) 

,  «X lot ,V*4A ( 120) 

COMMON  /SNIP/  4AS3.cINT,(14,c£*CE2«CE3,0MU.E0MU«e20MU.E30M0.t)F.eMM, 
.  NL  «FL, I A.E ( 1 20 ) 

COMMON  /CONST/  NCO. ECO. HI, OPR, RPU.UKAtfTr. RHO, NUM. MA(120) .CO.TA* 

H ( 120) .BtT A.HW(120) .TZ.ORAG.W.XO.T.XP.M.IT. 

OELTAS.TX.EST (120) ,NAR,MMAX(120) .TEST (120) , 

.  N  ( 120)  .PMALF 

COMMON  /In/  n* (120), 81(120) , VEL IN 
COMMON/OUT /NHR InT .NPLOT .ENO 

COMMON  /»AvEl/  2MA,Z4MA«EMAS.22WMA,22EMA,22oMA,E2MA2, 

1  ZwOOy (120) .ASlNPT (160.20) ,ACOSPT( 120.20) .CX6.SX6 

COMMON  /wave 2/  wO(20) ,K (20 ' «C(20) ,R«U (20) ,RW02 (20) ,RK (20) . 

.  R0(20) ,RWK(20) .M.HMS(20) 

COMMON  /INTER/  II.KTT(IO)  .01FK10) 

COMMON/ THAN5/ST  APT .RISE ,RAMH 
COMMON  /TtST/  VMA 

initialize  integral  sums 


154 

160 

161 

162 

163 

164 

165 

166 
167 

2 


1 


C  *  • 


MASS  *  0.0 
0 A  «  0.0 
IA  *  0.0 
CE  «  0.0 
CE2  =  0.0 
UMu  =0.0 

t0MU=0.0 
E20MU  «  0.0 
t3QMU  a  0.0 
«F  =  0. 0 
hmm  =  0.0 

2MA  =  0.0 
ZWMA  *  0.0 
fcMAS  =  0.0 


‘4V; 


'’**■*35* 


#0T 


ZZWMA  m 
ZWEMA  « 
Z2WMA  = 
t2MAZ  = 
VP4RT  = 
•  *  SET 


0.0 

0.0 

0.0 

0.0 


X ( 1 ) «SlN ( X(6)).X(2I »CUS(X(6> ) 
UH  FUNCTIONS  for  integrals  *  < 
UO  90  1  =  1.  -,1)M 
1PART ( I ) =t ( I ) *t ( I ) *M4 ( I > 

OPART ( I ) =t ( I ) *h a ( I ) 

ZWOT  =  0.0 


OAUX 
OAUX 
OAUX 
OAUX 
OAUX 
UAUX 
OAUX 
OAUX 
OAUX 
FunCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FunCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FUNCT 
FuNCT  46 
FUNCT  47 
FUNCT 
FUNCT 
FUNCT 
FUNCT 


3 

4 

5 

6 

7 

8 
4 

10 

11 

12 

13 

14 
16 
1 6 
17 
IN 
In 
20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

35 
37 

36 
34 

40 
4  1 

42 

43 

44 

45 


46 

49 

60 

61 


43 


U«OT  a  0.0 
0«0T  *  0.0 
li*0£  a  0.0 

XN£  a  *(J|.(N(ll*tX6-E<n*SX6) 

XNEO  a  X(J)*(E(1)*CX6*N(I) •SXO) 

UO  15  J  S  1.10 

ZrfQTT  «  -RrO(j)*AS1NPT(1.J> 

2«0T  *  ZaOT»2aCTT 

OVOTT  a  -Rr02Ij)*AC0SPTU.J> 

U*OT  a  D«tUT»OaOTT 

URDTT  a  -«MJ)*AS1NPT(I»J)»IAU)»C(J)«XnE) 

UHOT  *  OROT »ORGTT 

D*OE£  a  -RwKU)«ACOSPT<I.J> 

UVOE  a  OaOE.OaUtt 
IS  CONTINUE 

14 OOT(I)  a  /WUT«oaMP 
UVDE  *  OaOt»K‘*i«K 
URDT  *  OHOT*tiAAH 

U  a  X(l)*C>l6»X<2>*SK6.Z«0OT(l)*SAb 

VEL  *  VPAKT-X(  J)*t(I)-2aOOT(I>*Cx«» 

Z1(I)  a  MA(!)*Z«DOT (I) 

Z2(I>  *  MAI  I)  »Ul»DT*RAMP 
Z3( I)  a  EfI)*2c(I> 

Z«*(ll  a  E (  I ) *21  1 1) 

Z5(I)  *  rtA(I)«0*0«fOE 

26  ( 1 1  *  tlll'/iUI 

Z7|I)  «  *A(l)«vEL*U 

IF  I VEL.LE.O. >  00  TO  60 

IF  (81(11  .LE.O.U)  00  TO  50 

Ul(I)  a  OtL*Hl ( I)»(X<2)-XNEO-0R0T) 

00  TO  51 

50  ol(I)  a  0. 

51  CONTINUE 

l)2  <  I )  a  Ell) *0111) 

C1(I>  a  VtL*VtL*K  1) 

C2(  I  >  a  t ( t ) *Cl I f ) 
uO  TO  61 

60  01<I)  a  0. 
u2|l>  a  0. 

Cl (I>  3  o. 
t2<n  *  0. 

6 1  CONTINUE 

03 ( I )  a  «A ( I >  *u*0t* VEL 

04 1 1  >  a  Ell  )  *0J  1 1 

PIN  *  PI/2. 

05 ( 1 1  a  b(I)*(R«(I)-6(l)*TA/2«> 

66  06(1)  a  Ua ( 1 ) *t 1 1 >  *.5 

90  CONTINUE 
RrtOG=MFl'j*G3AVT  Y 
PlH  =  PI/2. 

H.PJ  a  K4R*PI 

Evacuate  intfgrals  using  trap  nt thoo 

1  a  1 
lNOEX  a  l 

91  CALC  THAFHAI  INuEX)  ,OIFF  1 1)  iMT  1 1 )  .TNASS) 
CALC  TRAP  I  'PAhI  (  InOEx)  .OIFFI  I)  .KTT  (1)  .OAl) 

calc  trap ici t t»utx> .oiffid  «ktt  m  «ceai 

CALC  Trap  (C2  •  li.uE  A)  .OIFF  ( I )  . M  T  ( 1 1  .CE2A) 
CALC  IRAKI  IPAhI  (  InUEX)  ,0  IFF  111  .MI  (1)  ,IAA> 


FUNCT  52 
FUNCT  53 
FUNCT  54 
FUNCT  55 
FUNCT  56 
FUNCT  5? 
FunCT  5tt 
FUNCT  59 
FUNCT  60 
FUNCT  61 
FUNCT  62 
FUNCT  63 
FUNCT  64 
FUNCT  65 
FUNCT  66 
FUNCT  6T 
FUNCT  68 
FUNCT  69 
FUNCT  70 
FUNCT  71 
FUNCT  7 2 
FUNCT  73 
FunCT  74 
FUNCT  75 
FUNCT  76 
FunCT  77 
FUNCT  78 
FunCT  79 
FUNCT  80 
FUNCT  51 
FUNCT  52 
FUNCT  A3 
FUNCT  84 
FUNCT  85 
FUNCT  86 
FUNCT  87 
FUNCT  A6 
FUNCT  A9 
FUNCT  90 
FUNCT  91 
FUNCT  92 
FUNCT  93 
FUNCT  94 
FunCT  95 
FUNCT  96 
FUNCT  9T 
FUNCT  90 
FunCT  59 
FUNCT 100 
FUNCT101 
FUNCT 102 
F  UNCT 103 
FunCT 1 04 
F  UNCT 1 05 
FUNCTION 
FUNCT 107 
FUNCTION 
FUNCT I  09 
FUNCT  UO 


CALL  TRAP (01 (INDEX! ,DIFF(I)«KTT(1) *DMUA) 

CALL  TRAP (02 (INDEX! .OIFF ( I ) .KTT ( 1 1 .EDMUA) 

CALL  TRAP(03(1NWLX) *OIFF ( l! .XT T ( 1 1 .E2DMUA) 

CALL  TRAP  ( 04  ( INDEX!  .OIFF ( 1 1 .KTT ( 1 >  .E30MUA) 

CALL  TRAP(l)5(  INOtXI  ,QIFF(I),KTT<1>  .BFA) 

CALL  TRAP (06 (INDEX) ,D1FF (1) .KTT ( 1 ! .8MMA) 

CALL  TRAP (71 (INDEX! .OIFF ( 1 ) «K TT ( l ) .ZHAA) 

CALL  TRAP (72 (INDEX! .OIFP ( 1) .KTT ( ll .ZwMAA) 

CALL  TRAP ( ?3 < INDEX) ,D1FF ( 1 ) .KTT ( I) .EMASA) 

CALL  TRAP ( 74  ( INDEX) ,DIFF<1) .KTT (I) . ZZwMAA) 

CALL  TRAP ( 7S( tNOEA) ,OIFF t l > ,KT 1(1 ) .ZwEMAA) 

CALL  TRAP ( 76 ( INDEX! .OIFF ( ll .KTT ( I ) , Z2WMAA) 

CALL  TRAP ( 77 ( INDEX! ,OIFF(I| ,KTT(1| . E2HAZA) 

93  CONTINUE 
HAsS  »  MASS  ♦  TNASS 

DA  «  QA  ♦  QAl 
IA  a  1A  .  I AA 
CE  *  CE  ♦  CEA 
CE2  «  CE2  .  CtcA 
UMj  «  OMU  .  OH(j  A 
EOMU  *  EOllI  ♦  EDMUA 
E20MU  *  E2')MU  ♦  EkUMUA 
E30MU  *  E30MU  ♦  E3DMUA 
BF  »  BF  *  jhugwBFm 
BN*)  *  BMM  ♦  RhdG»uMHA 
ZMa  *  ZMA*?MAA 
ZWMA  *  Z«MA*Z«*AA 
EMAS  *  EMAS*EMASA 
ZZWNA  «  ZZ4MA*ZZWMAA 
ZWENA  »  ZWEMA*2wEMAA 
Z2WMA  «  Z2.MA.22WHAA 
tZMAZ  «  E2  iA2*E2MAZA 

94  CONTINUE 
IF  (  I.GE.IIlOU  TO  92 
INDEX  »  iN'iEX.KTT  ( l )  -1 
1  *  1*1 
GO  TO  91 

92  CONTINUE 

•••••••  CALL  COMPUT  TO  FIND  THE  VALUE  OF  NL  AND  FL  USING 

THE  VALUES  OF  Th£  ABOVE  INTEGRALS 
CALL  COMPUT (X) 

I F ( NPR INT.LT.3)  GD  TO  111 
IF ( NPR INT  «E0. 3)  GD  TO  108 
IF  (NPR  I  NT  •  *-.  Q.**)  GO  TO  108 
WRIT£(6.n7»  (l*-ART<I),Isl,NUM! 

WRITE  (6,98)  (OPART(I),IaltNUM) 

WRITE (6.99)  ( C 1 ( I > • 1  *  1 .HUM) 

WHI TE (6, 1 00 )  (C2( 1 ) ,1=1 ,NDM) 

WRITE (6. 10?)  (Ul(l),l*l»NUM> 

WRITE  <6, 1 03)  (02(1) ,1  =  1 f NUM) 

WRITE  (6. I Oa )  (03(1)  ,1  =  1. NUM) 

WRITE  (6,1  OS)  (UA<1)  ,I  =  1»NUM> 

WRITE (6, 106)  IU6(I) .Ixl.NJM) 

WRI TE (6. 1 i 2)  (C6(l) ,I=1«NDM) 

WRITE (6. 11 3) <21 (1), I >1, NUM! 

WRITE (6.11a l Uc (ll ,I»1,NUM! 


FUNCTlll 
FUNCT 112 
FUNCT113 
FUNCT 114 
FUNCTllS 
FUNCT 116 
FUNCT 117 
FUNCT118 
FUNCT 119 
FUNCT 120 
FUNCT 121 
FUNCT 122 
FUNCT 123 
FUNCT 124 
FUNCT 125 
FUNCT 126 
FUNCT 127 
FUNCT 128 
FUNCT 129 
FUNCT130 
FuNCT 131 
FuNCT 132 
FUNCT133 
FUNCT 1 34 
FuNCT 135 
FUNCT 136 
FUNCT l 37 
FUNCT 138 
FUNCT 139 
FUNCT 140 
FUNCTlAl 
FUNCT1A2 
FUNCT 143 
FuNCT Iaa 
FUNCT 1 AS 
FUNCT >46 
FUNCTlAT 
FuNCT 148 
FUNCT 149 
FuNCT ISO 
FUNCT 151 
FUNCT 152 
FUNCT 153 
FUNCTISa 
FUNCT 155 
FUNCT 156 
FUNCT 1 57 
FUNCT 158 
FUNCT 159 
FUNCT 160 
FUNCT l6l 
FUNCT 16? 
FUNCT 163 
FUNCT  16*» 
FuNCT 165 
FUNCT 166 
FUNCTlhT 
FUNCT 168 
FUNCT 169 


WRITE <6.1 IS) (23(1). I *1 ,NUM) 

■RITE  (6.1  Is)  <2‘*(I)  .1*1. NUN) 

■RITE(6.1M)  (2S(  I)  .1*1, MUM) 

■RITE (6. 1 19) (26(1) .1=1, NUH) 

■RITE (6. 120) (2/ II) .1=1 .NOD 
■RITE(6.10/)KPI.Rm0g.PIh 
108  ■£ I TE  <6. 1 03)  MASS.CINT.0A.CE.CE2.CE3 
•RITE (6.121) [A 
121  FORMAT ( *  IA  •.ElO.A) 

■RITE(6.110)0M0.EL>Mu,E2fc)MU.E3UMU.»F,BMM 


FUNCT 170 
FUNCT1 71 
FUNCT 172 
FUNCT 173 
FUNCT 174 
FlINCT  1  7S 
FUNCT1  7b 
FUNCT 177 
FUNCT1 70 
FUNCT1 73 
FUNCT1R0 


96 

FORMAT (” 

CPAHT ( 1) ", 10(2X,E10.4> ) 

FUNCT181 

97 

FORMAT (” 

IPApI (I)”,lO(2X<£10.4) | 

Funct 182 

98 

FORMAT (” 

'JPAMI  (I)»,10(2X,E10.4>) 

FUNCT183 

99 

FORMAT <» 

Cl 

" *  10 (2X.E10.4) ) 

funct 1B4 

100 

FORMAT (" 

C2 

",  10 (2X.E1U.4) ) 

FUNCT185 

101 

FORMAT (“ 

C3 

" , 10 (2X.E10.4) ) 

FUNCT 106 

102 

FORMAT (” 

01 

",  10 (2X.E10.4) ) 

FUNCT187 

103 

FORMAT (" 

02 

"  ,  1 0 (2x  »E 1 0.4 ) ) 

FunCT 186 

104 

FORMAT (“ 

UJ 

",  10 (2X.E1U.4)  ) 

FUNCT 169 

106 

FORMAT (“ 

J4 

",10<2X,£i0.4| ) 

Funct i 90 

10b 

FORMAT (“ 

>)5 

" , 1 0 ( 2X ,E 10 .4) ) 

FUNCT 191 

112 

FORMAT (” 

06 

",  10 (2X.E10.4)  ) 

funct 192 

107 

FORMAT (» 

KPrtl 

".E1U,4."MH0G  ",E1 0.4,”  PMIM  M,E10.4) 

FUNCT 193 

109 

FORMAT (” 

mass 

".£10,4."  C1NT  ”,E10.4,»  QA  ",E10.4,”  CE  ”,£10.4, 

FUNCT 194 

•”Cc2  “,tlu.4.” 

CE3  ”,£10.4) 

FUNCT195 

110 

FORMAT <” 

JMU  " 

.E1U.4,”  EOMU  »,E1 0 .4,"  E20MU  ”.E10,4."  E30MU  ”. 

FUNCT196 

•El 0.4,”  bF  ”.£1u.4.”  BMM  ”,E10.4) 

113  F0RMAT(4M  71  ,10(2X,E10.4)) 

114  E0RMAT(4M  ?2  »10(2X,E10.4) ) 

115  FORMAT (4rt  23  .10 (2X.E10.4) ) 

116  F OrMAT (4m  24  .10(27, E10. 4)) 

118  F0RMAT(4M  25  .10(2*. E10. 4) ) 

113  FORMAT  (4r(  26  .  10  (2*,E10.4)  ) 

120  F  ORMAT  (4Fi  27  .10(27, E10.4)) 

117  FORMAT  (Sri  2MA  ,E10.4,bM  2aMA  ,E10.4,brf  EMAS  .E10.4, 

.  7h  22«ma  ,E10.4,7m  ZbEMa  .E10.4.7F)  ZZbma  .E10.4, 

.  7h  E2*A l  ,E10.4) 

111  CONTINUE 

mETURN 

ENO 

SUBROUTINE  COMPUT (X) 

DIMENSION  7(b) 

HEAL  KAR.KPI 
RE AC  NL.MASS.MA 

COMMON  /SHIP/  MASS«ClNT,QA,CE,CE2f CE J,L)MU»E0MU,E2OMU.E3DMU.BF ,BMM 
.  NC  ,FL  .  I  A  ,E  ( 1  20  ) 

COMMON  /CONST/  NCO.ECG.MI  ,Of'R,F(PO.G«AVTY,RMO,NUM,MA  (120)  .CD.TA, 

.  M  ( 120)  ,  HL  7  A  ,  FtW  (120)  ,7Z,l)PAG.».XU,T,XP,M,  IT. 

1  DECT AS. TX.EST (120) .KAR.MMAAI 120) .TEST (120) • 

,  n ( 120) .PmACF 

COMMON/OUT /NPM  1  NT  , NPLOT  , E NO 
COMMON  /TtOMS/  ll»T2»T3«T4,TS»Tb»T7»T8 
COMMON  /WAVE  1 /  ZMA .ZWMA.EMAS.ZZWMA, Z«EMA, Z?«MA,E?MAZ , 

1  ZbOOT (120) . AS 1NMT( 120.20) .ACUSPT ( 120,20) .CXb,SX6 

COMMON  /TEST/  OMA 


SEE  PAGES  7,6.  ANJ  9  OF  NOTES 


F  UnCT 197 
FUNCT 198 
FUNCT 199 
FUNCT  200 
FUNCT201 
FUNCT202 
FUNCT  20  3 
FUNCT  204 
FUNCT  205 
FUNCT206 
FUNCT207 

funct  208 

FUNCT209 
FUNCT  2 1 0 
COMPul  2 
COMPOT  3 
COMPUT  4 
COMPUT  b 
COMPUT  6 
.COMPul  7 
COMPUT  8 
COMPUT  9 
COMPUT 10 
COMPUI 1 1 
COMPUT 12 
COMPul 13 
COMPijT  1 4 
COMPUT 15 
compu  ns 
COMPUT 1 7 
COMPUT 18 
COMPUT 19 
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CX6  *  COS ( X (6) 1 

COMPUT20 

C0MPUT2I 

SX6  «  SIN(MO)) 

C0MPUT22 

PIN  a  Pl/2.0 

C0MPUT23 

KPf  a  KAR«PI 

COMPUT2** 

C0NS1  a  CX6 

C0MPUT25 

C0NS2  a  (NPI*RhO«PlH/TA)/CX6 

C0MPUT26 

CONSJ  a  CXN*SXo 

COMPUT27 

C0NS4  a  CX6*CXO 

C0MPUT2B 

TERN1  a  X<l>«CXb 

C0MPUT29 

TERN2  a  X<2)«S*0 

COMPUT30 

UV»|UM  a  <X<l)*CA6-<X(2>-Z*0OT<NUM) >*SX6)» 

COMPUT31 

4  <X< ll*SX6-X<3)*E<NUH>MXi2I-Z«00TTNuM>  >»CXfe> 

C0MPUT32 

ZMA  a  ZMA*X (3) *SX6 

CONPUT33 

C0MPUT34 

ZZWNA  a  ZZ4NA»A<3»*SX6 

COMPUT 35 

ZMHA  a  2*MA*CUNS1 

COMPUT  36 

ENAS  •  EMAS*CONSl 

COMPUT37 

UH(J  a  DMU*COnSa 

COMPUT 3B 

EOMU  *  Et)M()»CONSZ 

C0MPUT39 

CE  *  CE«CU*«PtvJ 

C0MPUT40 

CE2  «  CE2*cO«RRO 

COMPUT41 

t20«U  =  £20MU»C0NSJ 

COMPUT42 

E30MU  *  t30MU»CUNS3 

C0NPUT43 

2W£MA  a  2aEMA*CUNS4 

C0MPUT44 

Z2*NA  a  Z2*NA*CUNS4 

C0MPUT45 

20  T1  a  qa*x ( 3) *  <  T ERNI-TERN2 J 

C0MPUT46 
CONPUT  4  7 

T 1  a  T 1  ♦  ZZVMA  -  ENAS 

C0MPUT48 

T2  a  EOMU 

C0NPUT49 

T3  *  CE2 

COMPUTSO 

T4  a  MAtNUMlatOAUMIaUVNUH  ♦  E2HAZ  ♦  E30NU  - 

Z2RMA  •  BMH 

CONPUT SI 

NL  a  Tl  .  r2  *  T3  ♦  T4  •  dMH 

COMPUT52 

T5  a  HXSS»X<3»*(Tt»H2-TEft,»l> 

CONPUT53 

f5  »  T5  •  7*N A  -  /NA 

C0MPUT54 

T6  a  -OMU 

COMPUT55 

TT  a  -CE 

C0MPUT56 

T8  a  -MA(NUN)*UVNUH  -  E20NU  ♦  ZrEMA 

COMPUTS7 

BF  *  BF/CXl 

conputsb 

FLaT5«T6*T7*T8-bF 

CONPUT 59 
C0MPUT60 

1F»NPWINT.LT.3)G0  TO  30 

CONPUT61 
CONPUT 62 

25  CONTINUE 

C0NPUT63 

*RITET6*1U)NL*FL 

CONPUT 64 

10  FORMAT!"  NL  a  »*£12.6«"  FL  a  »,ti2.6) 

C0NPUT6S 

30  RETURN 

C0NPUT66 

end 

CONPUT67 

SUBROUTINE  INPUT 

INPUT  2 

••••••  DEFINITION  OF  INPUT  VARIABLES 

INPUT  3 

XA  a  INITIAL  TiM£ 

INPUT  4 

XE  *  FINAL  TIME 

INPUT  5 

MMJN  *  MINIMUM  STEP  SIZE 

INPUT  6 

UMAX  «  MAXIMUM  STEP  SIZE 

INPUT  7 

tPse  a  OELaT I Vt  ErrOR  CRITERION  USED  FUr  VALUES  OF  Y  6T  * 

INPUT  8 

EPS  *  ERROR  CRITERION  IN  KUT«tR 

INPUT  9 

A  a  ABSOLUTE  tRWOR  CRITERIA  USED  IN  KUTMER 

INPUT  10 

SPRINT  a  1  FINAL  PRINTOUT 

INPUT  11 

*  2  NaTRIX  INVERSE  MATRlX.F  COLUMN 

MATRIX* AND  KUTMER 

INPUT  12 

47 


.  -d.  *»•»'*•  v^i,^ .  -JnrMMH 
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RESULTS 

*  3  integral  values 

»  4  calculated  values-constant  for  given  inrut  values 

NPLOT  «  0  NU  PLOT 

a  1  PRINTER  PLOT 
End  «  number  of  runs 

N  *  MASS  OF  CRAFT 
V  *  WEIGHT  OF  CRAFT 

tz  =  thrust  component  in  l  direction 
tx  *  thrust  component  in  x  direction 
XECG  a  DISTANCE  From  CG  TO  CENTER  Of  PRESSURE  FOR  NORMAL  FORCE 
XP  *  MOMENT  arm  OF  PROPELLER  THRUST 

xo  *  distance  from  cg  to  center  uf  pressure  for  drag  force 

KAtll*  ADDED  mass  COEFFICIENT 

AN  ARRAY  GIvEN  THE  VALUE  (UR  mHICH  IS  READ  IN 
UM<II  a  -,EAM  At  FRtE  SURFACE  OR  AT  CHINE 

drag  *  FRICTION  drag 

K  a  WAVE  NUMBER 
RO  a  WAVE  HEIGHT 
NU  »  WAVE  SLOPE 
NU*t  =  NOMRER  OF  STATIONS 
8L  »  BOAT  LENGTH 
LAMBDA  a  WAVE  LENGTH 

RG  a  RADIUS  OF  GENERATION  IN  FEET 
T  *  PROPELLED  THRUST  IN  LHS 

Gamma  «  propeller  Thrust  anoLE  in  DEGREES 
OElTAS»STATIOn  SPACING  in  feet 

ECG  a  LONGITUDINAL  CENTER  OF  GRAVITY 
NCG  a  VERTICAL  LG 
BETA!  I)  a  DEAD  RISE 
NO ( I )  a  HEIGHT  OF  mean  BUTTOCN 
rHO  a  DENSITY  OF  baTER 
GHAVTY  a  GRAVITY  F1/SEC*«2 
DPR  a  DEGREES  PER  RADIAN 
rPU  a  RADIANS  per  OEGREE 

PI  *  3.14159  . 

t$T (II  a  STATION  POSITION 

START  *  START  TIME  OF  THE  RAMP  FUNCTION  FDR  SEA  WAVE 

RISE  *  DURATION  OF  THE  RISE  FROM  ZtRO  TO  ONE  OF  THE  RAMP 

•••••••1C  OPTIONS 

IC(l>  *1  USE  WAVE  Z  OISTANCE  IN  COMPUTING  LIFT  COMPONENT 
OF  NL  AND  FL 


REAL  IT  »A»M*MAf MMAX*NU»N*NCGtNO*MASS*NL* I  A* AAR 

Integer  end 


INRUT  13 
INPUT  14 
INPUT  15 
INPUT  16 
INPUT  17 
INPUT  IB 
INPUT  14 
INPUT  20 
INPUT  21 
INPUT  22 
INPUT  23 
INPUT  24 
INPUT  25 
INPUT  26 
INPUT  27 
INPUT  2B 
INPUT  29 
Input  30 
INPUT  31 
INPUT  32 
INPUT  33 
INPUT  34 
INPUT  35 
INPUT  36 
INPUT  37 
INPUT  38 
INPUT  39 
INPUT  40 
INPUT  41 
INPUT  42 
INPUT  43 
INPUT  44 
INPUT  45 
INPUT  46 

Input  47 
input  48 

INPUT  49 
INPUT  50 
INPUT  51 
INPUT  52 

input  53 

INPUT  5A 
INPUT  55 
INPUT  56 

Input  57 

INPUT  58 
INPUT  59 
INPUT  60 
INPUT  61 
INPUT  62 


INPUT  63 

COMMON  /CONST/  NCG,£CG«PI»DPR*RPO»GRAVTY,RHO*NUM,MA(1201.CO»TA,  INPUT  64 
«  B(12u> (BETA. Hwf 120) .TZ.ORAG»W.XD.T.XPfM«IT*  INPUT  65 
1  DELTAS, TX,EST(120) .«AR,MM«X(120> .TEST(120) ,  INPUT  66 
.  NU20»,PHALF  INPUT  67 
COMMON  /SHIP/  MASS«CINT *0 a*CE ,CE2 ,CE3,OMU,EOMU,E2DMU,E3DMU»0F ,  0MM, INPUT  6fl 
*  NL  *FL « I A  «E  ( 1  <?0 )  INPUT  69 
COMMON  /IN/  WM ( 120) , y 1 (120) fVELIN  INPUT  70 
COMMON  /IN?/  NC(12U| ,XA,XE.HRAXfHMlN,A(61 *EPSt«6>  INPUT  71 
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COMMON/OUT /NPRlNT.NPLOTlE NO 
common  /ACCEL/  XACCL,8WACl«CGACL«6L 

NAmELIST/HSP/A.MPrInT .NPLUT .END.M.BL,  TZ.TX.XECG.XP.XD, 
t  DRAG, RG,T .GAMMA, ECG,NCG,KAR,NUM,BETA,EST 

.  ,XA,X£,HNIM»HMAX,ERS»9£LIN 

OAT*  A  /. 01, .0001,. 00081, .U.0001*. 00001/ 

OAT*  NPRINT.NPLOT, END/1, 1.1/ 

DATA  M.HL , TZ»  TX.XtCG.XP, XO.URAG.RG, T , GAMMA, 

1  tCG.NCG.i^R  /lb. ,3.  73, o*0. 0, .962b«2*0.0« 

.  2.325,0.0,1.0/ 

DATA  NUM,8£TA,tST  /77.20.0, 

4  0. 0000,. 03Ub,. 0b250,. 09375,. 12500,. 15625,. 18750,. 21875, 

.  .25000,. 28115.. 31250,. 3*376,. 3 7500,. 40625,. 43750,. 46875, 

.  . 50000. .531CS,. 56250, .39375, .62500, .65625,. 6875*. 71875. 

i  .  75000*.  7dU5,. 81250,  .84375,. 8  7500,.  <70625,  .93750*. 968  75, 1.000, 

.  1. 06250.1. 12500. 1.187S0, 1.26000,1.3125,1. 37500, 1.4375, 

.  1.500, 1.5625, 1.625, 1.6675, 1.75. 1.8125,1. *75, 1.9375, 2.0. 

.  2. 0625, 2. 125, 2. 1875, 2. 25, 2. 31 25, 2. 375, 2. *375, 2. 5, 2. 5625, 2. 625, 

,  2.6875,2.73.2.6125,2.8750,2.9375,3.0,3.0625,3.125,3.1875, 

.  3.2500.3.3125  ,3. 375,3. *373,3.5, 3. 5626, 3. 62S, 3. 6875, 3. 75  / 

DATA  XA. XE, HM1N, HMAX.EPS  /0. 0,20.0, .025,.!, .15/ 

DATA  VELIN  /19.62/ 

•  •••••••  BfcAU  IN  AND  aRlIE  OUT  KUTMER  PARAMETERS  ANO  PROGRAM 

OPTIONS 

HEA0(5.M5P» 

■RITE (6»MSP) 

00  10  1=1, n 
10  EPsEdl  *  EPS 


,•••••••  SET  UR  CONSTANTS 

PI  *  3.1M592653S69 
GRAVTY*32. 18 
OPR*57.2957795l30tt 
RRD“. 01 7453292519 
IF  (EST(NUM).LT.8LI  STOP  3 

COMPUTE  NO  ANO  RM  ARRAYS 

THIS  IS  FOR  SPtClAL  BOA  FORM  ONLY.  CHANGE  PROGRAM 
THRU  STATEMENT  32  FOR  NE6  BU*  SHAPE 

00  32  1  =  1,  NIJM 

IF (EST < I) .Gfc.O. 75)  GO  TO  30 

nO( 1 1 »-0.*6875« ( 1 • 0-SQRT (t5T(l)/0.375-(ESTCll/0.7Sj»*2.0l 1 
BM( 1 1 «.375«S0RT (1.0. (EST 111/. 76-1. 1 ••2.0> 

GO  TO  32 
30  NO (11=0.0 

oM( I )  »  0.375 
32  CONTINUE 

.•••••COMPUTE  CONSTANTS  ANO  INITIALIZE  ARRAYS 
M*R/GHAVTY 
RHQ=1 .99 
IT«M#RG*Ko 
PHALF  »  (PI/2.) *rhO 

aETA  *  8ETa*RPu 
CO  *  COS(dETA) 


INPUT  72 
INPUT  T3 
INPUT  7* 
INPUT  75 
INPUT  76 
INPUT  77 
INPUT  78 
INPUT  79 
INPUT  80 
INPUT  81 
INPUT  82 
INPUT  83 
INPUT  8* 
INPUT  85 
INPUT  86 
INPUT  87 
INPUT  88 
INPUT  89 
INPUT  90 
INPUT  91 
INPUT  92 
INPUT  93 
INPUT  94 
INPUT  95 
INPUT  9b 
INPUT  97 
INPUT  98 
INPUT  99 
INPUT  100 
INPUT  101 

Input 102 
INPUT  103 
INPUT  104 

Input  105 

INPUT  106 
INPUT 107 
INPUT l 08 
INPUT  109 
INPUT  1 10 
Input  ill 
Input 112 
INPUT  113 
INPUT  11* 
INP.iTl  15 
INPUT  116 
INPUT  1 17 
INPUT 116 
INPUT 119 
INPUT  120 
INPUT 121 
INPUT  1 22 
INPUT  123 

Input  12* 
input  125 
INPUT  12b 

Input  127 
INPUT  128 

Input  129 
Input  130 
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OOUOOUDUUO  u  u  uo  u  u  u  uu 


TA  >  TAN(BET4I 
UO  60  1*1. NUN 

E(I>  *  ECG-EST(l) 

N(  I)  •  NCG.NOUl 

MMAX(I)  *  KAR«PHALF*8M(ll*BM(ll 
TEST ( I >  *  (2.»e«(I»«TA)/Pi 
6«  CONTINUE 
ENO*EnO«1 
return 
ENO 

subroutine  ploter  (fx.xa.hhax.ib.ipti 
input: 

FX  A  TWO  DIMENSIONAL  ARRAY  CONTAINING  PITCH  AND 

HEAVE  VALUES  AT  EACH  TINE  STEP 

XA  INITIAL  TIME 

HNAX  TIME  INTERVAL*  MT1mE*HHAX  ■  INTERVAL  BETWEEN 

f *  VALUES 

Id  NUMBER  OF  FX  VALUES 


REAL  IT.K.LAMdUA.H.MA.MMAX.N.NCG 
INTEGER  ENT 

OImENSION  F X  (2,2000)  fFMIN  (2)  «f  MAXI2)  ,NVAR(2) 

COMMON  /CONST/  NCG.ECG«MI*DHR.RPU.GRAVTY.RHJ,nUM.MAc120).CD.TA, 
.  d 1 120)  .BETA. n«i(  1201  ♦  T2, DRAG. i*,XD.T.XP.M.  IT, 

1  DELTAS. TA.EST  (120) .KAR.MMAX ( 1201 .TEST (120) , 

.  N(120),PHALF 

COmMON/OUT/NPR I NT . NPLOT  «£NU 

••••••••  SET  UP  VALUES  FOR  PLOT  AND  CREATE  PLOT 

Nf  uN*2 

••*«*•••  SE I  UP  MIN  ANO  MAX  LIMITS  FOR  PLOT 

FMIN(i>zFX(l.i> 

FMIN(2»*FX(2.H 

FMAX(1)«FX(1.1> 

FMAX(2)«FX(2*II 

00  200  1*1. Id 

IF(FX(1,I) .LT.FMIN(I) JFMINU JafXfl.I) 

IF  (FX  <  I  ,  1 1  .GT.FMAX(lMFMAX<H*fX(i»I> 

1F(FX(2.I> .LT.fMlN(2l )FNlN<2)*Fx(2,I) 

If(FX(2,I).GT.f max(2))FMAX(2>»Fx(2.I> 

200  CONTINUE 

POO  CONTINUE 

NVAR ( 1 1 *10H  HEAVE 
NVAR (2)  =  1UH  PI  I CH 
Nl*2 
AO«XA 

UELX  *  HMAX 

IF  (NPLOT. E'l.nCALL  PLOT2(FX.FhIN»FMAX,NVAR,NFUN,N1,1B,XO.DELX> 
RETURN 

end 

SUBROUTINE  TRAP (F.UX.NPTS. ANSI 
INPUT* 


INPUT  131 
INPUT 132 
INPUT 133 

Input 134 

INPUT l 3S 

INPUT  136 

INPUT 137 

INPUT 13d 

INPUT 139 

INPUT 140 

PLOTEk  2 

PLOTER  3 

PLOTER  4 

PLOTER  5 

PLOTER  6 

PLOTER  7 

PLOTER  d 

PLOTEP  9 

PL0TEP10 

PL0TER11 

PLOTER 12 

PL0TER13 

PL0TER14 

PL0TER1S 

PL0TER16 

PL0TER17 

PL0TER18 

PL0TER19 

PL0TEP20 

PL0TER21 

PL0TER22 

PLOTER23 

PL0TER24 

PL0TER2S 

PL0TER26 

PL0TER27 

PL0TER2d 

PL0TER29 

PLOTER 30 

PL0TEP31 

PL0TER32 

PL0TER33 

PL0TER34 

PL0TER3S 

PL0TER36 

PL0TER37 

PL0TER3S 

PL0TER39 

PLOTER40 

PL0TER41 

PL0TER42 

PL0TER43 

RL0TERA4 

PL0TER45 

PL0TER46 

PL0TER47 

TRAP  2 

TPAP  3 

TRAP  4 


uuuuuu  uuuuuu 


F  ARRAY  Of  FUNCTIONAL  VALUES  Of  THE  INTEGRA*® 

OX  Tnt  X  INTERVAL  MtTwEtN  VALUES 

NPTS  THE  NUMBER  Of  VALUES  GIVEN 

OUTPUT  t 

ANS  The  value  Of  The  INTEGRAL 

OIHENSIOn  f(NPTS) 

ANSaO.O 

IF(NPTS.LT.2)GU  TO  99V 
1)0  1  la  1.  NPTS 
1  ANSaANS*MI) 

ANS«OXa(ANS-O.S«(f  <U»f  (NfTS>  )) 

999  CONTINUE 
RETURN 
END 

FUNCTION  RnP(T.START.RISE) 

•  •  •  •  •  •  •  THIS  f UNCTION  IS  USEO  TO  GRADUALLY  IHPLIHENT  THE  HAVE 


T  CURRENT  TINE 

START  TIHE  TO  START  RAMP  f huh  0.0  TO  1.0 

RISE  THE  LENGTH  Of  THt  RISE  fROM  0.0  TO  1.0 

Hao.o 

If (T.LT.START)GO  TO  99 
If  (RISE.K'J.O.O)uO  TO  80 
TOPaT-STAHT 

H*1 . 0 

If  (TOH.LT.RlSUH.TOP/RISE 
GO  TO  99 

SO  Hal. 

if (T.EO.STsRT)HaO.S 
99  WHPan 
RETURN 
ENq 

SUBROUTINE  SEAaAY 
REAL  K  . 

COMHON  /WAVE//  w0(20) ,K(20) «C(/0) *RwO(/Q) *RWO/(20) *RK(20) » 

*  R0(20>  .K.K  (20)  ,h,pmS(20) 

COMHON  /CONST/  NCG«ECG.PI *0PR*RPU*6RAVT Y*RHJ,NUM.MA  (120) .CO.TA. 

•  8(120) «bET4.hw( 120) * T2*DRAG*w«XDt T. XP.M, IT* 

1  0ELTAS.TX«EST(12u> .KAR.mMAXC 12U) .TEST  1120) t 

.  N(120),PH4LF 

READ  IS. 00)  H 
.RITE  (6. VI)  H 
80  FORMAT  (8F10.A) 

90  FORMAT  (2*m  SIGNIFICANT  aAVE  HEIGHT  a  ,E10.A,2X.5H  FEET,//) 
HROOT  *  SORT (h) 

■N  a  2.276/HRuCT 
aU( 1 )  «  ,/9SA**n 
■0(2)  a  vN 
■0(3)  «  1.1S3**N 
■0(A)  a  1.403b*wN 
■0(5)  a  1.60 1 9*WN 
■  0(6)  a  1 »  796 J*  »N 
■0<T)  a  2.0Q36««N 
■0(8)  «  2.19A1*4N 
■0(9)  a  2.7<»IV»WN 
■0(10)  a  2.6l2*aN 
R0(1 )  a  0.  I  3b**«H 


TRAP  S 
TRAP  6 
TRAP  T 
TRAP  8 
TRAP  9 
TRAP  10 
TRAP  1 1 
TRAP  *  12 
TRAP  13 
TRAP  1A 
TRAP  IS 
TRAP  |6 
TRAP  17 
TRAP  18 
Trap  19 
RHP  2 
RHP  3 
RHP  4 
MHP  S 
RHP  6 
RHP  7 
RHP  8 
RHP  9 
RHP  10 
RHP  11 
RHP  12 
RHP  13 
RHP  1A 
RHP  IS 
RHP  16 
RHP  17 
RHP  18 
RHP  19 
RHP  20 
SEAWAY  2 
SEAWAY  3 
SEAaAY  4 
SEAWAY  S 
SEAWAT  6 
SEAWAY  7 
SEAWAY  8 
SEAWAY  V 
SEAWAY I  0 
SEAWAY! 1 
SEAWAY) 2 
SE  AW AY  1 3 
SEAWAY14 
SEAWAV15 
SEAWAY  16 
SEAWAY17 
SEAWAY  1 8 
SEAWAY19 
SEAWAY20 
SEAWAY!  1 
SEAWAY?/ 
SEAWAY  23 
SEAWAY/4 
SEA.AY2S 
SE AWAY/b 


i 


K0|2)  a  0.1861** 

«0|3)  a  0.1657*ri 
R0<4)  a  0.13U2«H 
ROjS)  a  0.099<*»rl 
*0(61  a  0.u77l»H 
HO  ft)  a  O.O604«H 
*0(8)  a  0..14»2»rt 
*0(9)  a  U.039U»rl 
RO(IO)  a  U.062b*H 
PHSIH  «  .005 
♦’MS*-*)  *  2.41 
P«S<3»  a  5.20 
HHS(41  a  4.00 
P«SIS»  a  0.00 
»*MS<6)  a  1.27 
HNS!/)  a  3.11 
**«S<8I  a  2.92 
*'HS(9l  a  3.55 
t’HSUO)  a  0.70 
UO  50  Jal.io 

Mj)  a  W0(  jl*»0  ( Jl  /GRAVTY 
C(J)  a  GRAVTV/aO<J) 

RHOU)  a  i?0(J)  «aO(  J) 
RM02IJI  a  300(.>)»a0(j) 
a  W«0  (  J)  *K  <  j) 

HK(J)  a 

50  CONTINUE 
RETURN 
tNO 

SUBROUTINE  TABLE 

CUMNON/S I NF /RO I n T I 1 0 0 0 ) 

U*  a  ,003l4l592t»54 
*  a  1.570746347 
UO  100  Jal.501 

POINT  f J)  a  S1N(X) 


*  a  1002-J 

POINT (Ki  * 

X  a  x«0X 

100  CONTINUE 
RETURN 
ENO 


*  -HUINI ( J) 


SEAHAY27 

SEAHAY28 

SEAHAY29 

SEAMAY30 

SEAHAY31 

SEAHAV32 

SEAMAY33 

SEAHAY34 

SEANAY3S 

SEANAY36 

SEAHAY37 

SEAHAY38 

SEAMAY39 

SEAhA Y4Q 

seaRayai 
SEAMA Y42 
SEAHAY43 
SEAHAY44 
SEAWA Y4S 
SE4mAY46 
SEAUAY47 
SE AaA Y48 
SEAWA Y49 

seawayso 
SEAHA YSl 
SEAHAY52 
SEAaA Y53 
SEAHAY54 
SEAHAY55 

table  2 
table  3 

TABLE  4 

table  5 
TABLE  6 
TABLE  7 

table  8 
table  9 
table  10 
table  11 
table  12 

TABLE  13 


DTNSRDC  ISSUES  THREE  TYPES  OF  REPORTS 

1.  DTNSRDC  REPORTS.  A  FORMAL  SERIES,  CONTAIN  INFORMATION  OF  PERMANENT  TECH 
NICAL  VALUE.  THEY  CARRY  A  CONSECUTIVE  NUMERICAL  IDENTIFICATION  REGARDLESS  OF 
THEIR  CLASSIFICATION  OR  THE  ORIGINATING  DEPARTMENT. 

2.  DEPARTMENTAL  REPORTS,  A  SEMIFORMAL  SERIES.  CONTAIN  INFORMATION  OF  A  PRELIM 
INARY,  TEMPORARY.  OR  PROPRIETARY  NATURE  OR  OF  LIMITED  INTEREST  OR  SIGNIFICANCE. 
THEY  CARRY  A  DEPARTMENTAL  ALPHANUMERICAL  IDENTIFICATION. 

3.  TECHNICAL  MEMORANDA.  AN  INFORMAL  SERIES,  CONTAIN  TECHNICAL  DOCUMENTATION 
OF  LIMITED  USE  AND  INTEREST.  THEY  ARE  PRIMARILY  WORKING  PAPERS  INTENDED  FOR  IN 
TERNAL  USE.  THEY  CARRY  AN  IDENTIFYING  NUMBER  WHICH  INDICATES  THEIR  TYPE  AND  THE 
NUMERICAL  CODE  OF  THE  ORIGINATING  DEPARTMENT.  ANY  DISTRIBUTION  OUTSIDE  DTNSRDC 
MUST  BE  APPROVED  BY  THE  HEAD  OF  THE  ORIGINATING  DEPARTMENT  ON  A  CASE-BY-CASE 
BASIS. 


